Actual source code: iterativ.c

  1: /*
  2:    This file contains some simple default routines.
  3:    These routines should be SHORT, since they will be included in every
  4:    executable image that uses the iterative routines (note that, through
  5:    the registry system, we provide a way to load only the truly necessary
  6:    files)
  7:  */
  8: #include <petsc/private/kspimpl.h>
  9: #include <petscdmshell.h>
 10: #include <petscdraw.h>

 12: /*@
 13:    KSPGetResidualNorm - Gets the last (possibly approximate and/or preconditioned) residual norm that has been computed.

 15:    Not Collective

 17:    Input Parameters:
 18: .  ksp - the iterative context

 20:    Output Parameters:
 21: .  rnorm - residual norm

 23:    Level: intermediate

 25:    Note:
 26:    For some methods, such as `KSPGMRES`, the norm is not computed directly from the residual.

 28:    The type of norm used by the method can be controlled with `KSPSetNormType()`

 30: .seealso: [](chapter_ksp), `KSP`, `KSPSetNormType()`, `KSPBuildResidual()`, `KSPNormType`
 31: @*/
 32: PetscErrorCode KSPGetResidualNorm(KSP ksp, PetscReal *rnorm)
 33: {
 36:   *rnorm = ksp->rnorm;
 37:   return 0;
 38: }

 40: /*@
 41:    KSPGetIterationNumber - Gets the current iteration number; if the `KSPSolve()` is complete, returns the number of iterations used.

 43:    Not Collective

 45:    Input Parameters:
 46: .  ksp - the iterative context

 48:    Output Parameters:
 49: .  its - number of iterations

 51:    Level: intermediate

 53:    Note:
 54:    During the ith iteration this returns i-1

 56: .seealso: [](chapter_ksp), `KSP`, `KSPGetResidualNorm()`, `KSPBuildResidual()`, `KSPGetResidualNorm()`, `KSPGetTotalIterations()`
 57: @*/
 58: PetscErrorCode KSPGetIterationNumber(KSP ksp, PetscInt *its)
 59: {
 62:   *its = ksp->its;
 63:   return 0;
 64: }

 66: /*@
 67:    KSPGetTotalIterations - Gets the total number of iterations this `KSP` object has performed since was created, counted over all linear solves

 69:    Not Collective

 71:    Input Parameters:
 72: .  ksp - the iterative context

 74:    Output Parameters:
 75: .  its - total number of iterations

 77:    Level: intermediate

 79:    Note:
 80:     Use `KSPGetIterationNumber()` to get the count for the most recent solve only
 81:    If this is called within a linear solve (such as in a `KSPMonitor` routine) then it does not include iterations within that current solve

 83: .seealso: [](chapter_ksp), `KSP`, `KSPBuildResidual()`, `KSPGetResidualNorm()`, `KSPGetIterationNumber()`
 84: @*/
 85: PetscErrorCode KSPGetTotalIterations(KSP ksp, PetscInt *its)
 86: {
 89:   *its = ksp->totalits;
 90:   return 0;
 91: }

 93: /*@C
 94:   KSPMonitorResidual - Print the (possibly preconditioned) residual norm at each iteration of an iterative solver.

 96:   Collective

 98:   Input Parameters:
 99: + ksp   - iterative context
100: . n     - iteration number
101: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
102: - vf    - The viewer context

104:   Options Database Key:
105: . -ksp_monitor - Activates KSPMonitorResidual()

107:   Level: intermediate

109:   Note:
110:   For some methods, such as `KSPGMRES`, the norm is not computed directly from the residual.

112:   The type of norm used by the method can be controlled with `KSPSetNormType()`

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

117: .seealso: [](chapter_ksp), `KSPMonitorSet()`, `KSPMonitorTrueResidual()`, `KSPMonitorResidualDraw()`, `KSPMonitorResidualDrawLG()`,
118:           `KSPMonitorResidualRange()`, `KSPMonitorTrueResidualDraw()`, `KSPMonitorTrueResidualDrawLG()`, `KSPMonitorTrueResidualMax()`,
119:           `KSPMonitorSingularValue()`, `KSPMonitorSolutionDrawLG()`, `KSPMonitorSolutionDraw()`, `KSPMonitorSolution()`,
120:           `KSPMonitorErrorDrawLG()`, `KSPMonitorErrorDraw()`, KSPMonitorError()`
121: @*/
122: PetscErrorCode KSPMonitorResidual(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
123: {
124:   PetscViewer       viewer = vf->viewer;
125:   PetscViewerFormat format = vf->format;
126:   PetscInt          tablevel;
127:   const char       *prefix;

130:   PetscObjectGetTabLevel((PetscObject)ksp, &tablevel);
131:   PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix);
132:   PetscViewerPushFormat(viewer, format);
133:   PetscViewerASCIIAddTab(viewer, tablevel);
134:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
135:   PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " KSP Residual norm %14.12e \n", n, (double)rnorm);
136:   PetscViewerASCIISubtractTab(viewer, tablevel);
137:   PetscViewerPopFormat(viewer);
138:   return 0;
139: }

141: /*@C
142:   KSPMonitorResidualDraw - Plots the (possibly preconditioned) residual at each iteration of an iterative solver.

144:   Collective

146:   Input Parameters:
147: + ksp   - iterative context
148: . n     - iteration number
149: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
150: - vf    - The viewer context

152:   Options Database Key:
153: . -ksp_monitor draw - Activates `KSPMonitorResidualDraw()`

155:   Level: intermediate

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

161: .seealso: [](chapter_ksp), `KSPMonitorSet()`, `KSPMonitorTrueResidual()`, , `KSPMonitorResidual()`, `KSPMonitorResidualDrawLG()`
162: @*/
163: PetscErrorCode KSPMonitorResidualDraw(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
164: {
165:   PetscViewer       viewer = vf->viewer;
166:   PetscViewerFormat format = vf->format;
167:   Vec               r;

170:   PetscViewerPushFormat(viewer, format);
171:   KSPBuildResidual(ksp, NULL, NULL, &r);
172:   PetscObjectSetName((PetscObject)r, "Residual");
173:   PetscObjectCompose((PetscObject)r, "__Vec_bc_zero__", (PetscObject)ksp);
174:   VecView(r, viewer);
175:   PetscObjectCompose((PetscObject)r, "__Vec_bc_zero__", NULL);
176:   VecDestroy(&r);
177:   PetscViewerPopFormat(viewer);
178:   return 0;
179: }

181: /*@C
182:   KSPMonitorResidualDrawLG - Plots the (possibly preconditioned) residual norm at each iteration of an iterative solver.

184:   Collective

186:   Input Parameters:
187: + ksp   - iterative context
188: . n     - iteration number
189: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
190: - vf    - The viewer context

192:   Options Database Key:
193: . -ksp_monitor draw::draw_lg - Activates `KSPMonitorResidualDrawLG()`

195:   Level: intermediate

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

201:   Use `KSPMonitorResidualDrawLGCreate()` to create the contex used with this monitor

203: .seealso: [](chapter_ksp), `PETSCVIEWERDRAW`, `KSPMonitorSet()`, `KSPMonitorTrueResidual()`, `KSPMonitorResidualDraw()`, `KSPMonitorResidual()`
204: @*/
205: PetscErrorCode KSPMonitorResidualDrawLG(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
206: {
207:   PetscViewer        viewer = vf->viewer;
208:   PetscViewerFormat  format = vf->format;
209:   PetscDrawLG        lg     = vf->lg;
210:   KSPConvergedReason reason;
211:   PetscReal          x, y;

215:   PetscViewerPushFormat(viewer, format);
216:   if (!n) PetscDrawLGReset(lg);
217:   x = (PetscReal)n;
218:   if (rnorm > 0.0) y = PetscLog10Real(rnorm);
219:   else y = -15.0;
220:   PetscDrawLGAddPoint(lg, &x, &y);
221:   KSPGetConvergedReason(ksp, &reason);
222:   if (n <= 20 || !(n % 5) || reason) {
223:     PetscDrawLGDraw(lg);
224:     PetscDrawLGSave(lg);
225:   }
226:   PetscViewerPopFormat(viewer);
227:   return 0;
228: }

230: /*@C
231:   KSPMonitorResidualDrawLGCreate - Creates the context for the (possibly preconditioned) residual norm monitor `KSPMonitorResidualDrawLG()`

233:   Collective

235:   Input Parameters:
236: + viewer - The `PetscViewer` of type `PETSCVIEWERDRAW`
237: . format - The viewer format
238: - ctx    - An optional user context

240:   Output Parameter:
241: . vf    - The viewer context

243:   Level: intermediate

245: .seealso: [](chapter_ksp), `KSP`, `PETSCVIEWERDRAW`, `KSPMonitorSet()`, `KSPMonitorTrueResidual()`, `KSPMonitorResidualDrawLG()`
246: @*/
247: PetscErrorCode KSPMonitorResidualDrawLGCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
248: {
249:   PetscViewerAndFormatCreate(viewer, format, vf);
250:   (*vf)->data = ctx;
251:   KSPMonitorLGCreate(PetscObjectComm((PetscObject)viewer), NULL, NULL, "Log Residual Norm", 1, NULL, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &(*vf)->lg);
252:   return 0;
253: }

255: /*
256:   This is the same as KSPMonitorResidual() except it prints fewer digits of the residual as the residual gets smaller.
257:   This is because the later digits are meaningless and are often different on different machines; by using this routine different
258:   machines will usually generate the same output.

260:   Deprecated: Intentionally has no manual page
261: */
262: PetscErrorCode KSPMonitorResidualShort(KSP ksp, PetscInt its, PetscReal fnorm, PetscViewerAndFormat *vf)
263: {
264:   PetscViewer       viewer = vf->viewer;
265:   PetscViewerFormat format = vf->format;
266:   PetscInt          tablevel;
267:   const char       *prefix;

270:   PetscObjectGetTabLevel((PetscObject)ksp, &tablevel);
271:   PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix);
272:   PetscViewerPushFormat(viewer, format);
273:   PetscViewerASCIIAddTab(viewer, tablevel);
274:   if (its == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
275:   if (fnorm > 1.e-9) PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " KSP Residual norm %g \n", its, (double)fnorm);
276:   else if (fnorm > 1.e-11) PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " KSP Residual norm %5.3e \n", its, (double)fnorm);
277:   else PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " KSP Residual norm < 1.e-11\n", its);
278:   PetscViewerASCIISubtractTab(viewer, tablevel);
279:   PetscViewerPopFormat(viewer);
280:   return 0;
281: }

283: PetscErrorCode KSPMonitorRange_Private(KSP ksp, PetscInt it, PetscReal *per)
284: {
285:   Vec                resid;
286:   const PetscScalar *r;
287:   PetscReal          rmax, pwork;
288:   PetscInt           i, n, N;

290:   KSPBuildResidual(ksp, NULL, NULL, &resid);
291:   VecNorm(resid, NORM_INFINITY, &rmax);
292:   VecGetLocalSize(resid, &n);
293:   VecGetSize(resid, &N);
294:   VecGetArrayRead(resid, &r);
295:   pwork = 0.0;
296:   for (i = 0; i < n; ++i) pwork += (PetscAbsScalar(r[i]) > .20 * rmax);
297:   VecRestoreArrayRead(resid, &r);
298:   VecDestroy(&resid);
299:   MPIU_Allreduce(&pwork, per, 1, MPIU_REAL, MPIU_SUM, PetscObjectComm((PetscObject)ksp));
300:   *per = *per / N;
301:   return 0;
302: }

304: /*@C
305:   KSPMonitorResidualRange - Prints the percentage of residual elements that are more then 10 percent of the maximum value.

307:   Collective

309:   Input Parameters:
310: + ksp   - iterative context
311: . it    - iteration number
312: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
313: - vf    - The viewer context

315:   Options Database Key:
316: . -ksp_monitor_range - Activates `KSPMonitorResidualRange()`

318:   Level: intermediate

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

324:  .seealso: [](chapter_ksp), `KSP`, `KSPMonitorSet()`, `KSPMonitorResidual()`
325: @*/
326: PetscErrorCode KSPMonitorResidualRange(KSP ksp, PetscInt it, PetscReal rnorm, PetscViewerAndFormat *vf)
327: {
328:   static PetscReal  prev;
329:   PetscViewer       viewer = vf->viewer;
330:   PetscViewerFormat format = vf->format;
331:   PetscInt          tablevel;
332:   const char       *prefix;
333:   PetscReal         perc, rel;

336:   PetscObjectGetTabLevel((PetscObject)ksp, &tablevel);
337:   PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix);
338:   PetscViewerPushFormat(viewer, format);
339:   PetscViewerASCIIAddTab(viewer, tablevel);
340:   if (!it) prev = rnorm;
341:   if (it == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
342:   KSPMonitorRange_Private(ksp, it, &perc);
343:   rel  = (prev - rnorm) / prev;
344:   prev = rnorm;
345:   PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " KSP 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));
346:   PetscViewerASCIISubtractTab(viewer, tablevel);
347:   PetscViewerPopFormat(viewer);
348:   return 0;
349: }

351: /*@C
352:   KSPMonitorTrueResidual - Prints the true residual norm, as well as the (possibly preconditioned) approximate residual norm, at each iteration of an iterative solver.

354:   Collective

356:   Input Parameters:
357: + ksp   - iterative context
358: . n     - iteration number
359: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
360: - vf    - The viewer context

362:   Options Database Key:
363: . -ksp_monitor_true_residual - Activates `KSPMonitorTrueResidual()`

365:   Level: intermediate

367:   Notes:
368:   When using right preconditioning, these values are equivalent.

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

373: .seealso: [](chapter_ksp), `KSP`, `KSPMonitorSet()`, `KSPMonitorResidual()`, `KSPMonitorTrueResidualMaxNorm()`
374: @*/
375: PetscErrorCode KSPMonitorTrueResidual(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
376: {
377:   PetscViewer       viewer = vf->viewer;
378:   PetscViewerFormat format = vf->format;
379:   Vec               r;
380:   PetscReal         truenorm, bnorm;
381:   char              normtype[256];
382:   PetscInt          tablevel;
383:   const char       *prefix;

386:   PetscObjectGetTabLevel((PetscObject)ksp, &tablevel);
387:   PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix);
388:   PetscStrncpy(normtype, KSPNormTypes[ksp->normtype], sizeof(normtype));
389:   PetscStrtolower(normtype);
390:   KSPBuildResidual(ksp, NULL, NULL, &r);
391:   VecNorm(r, NORM_2, &truenorm);
392:   VecNorm(ksp->vec_rhs, NORM_2, &bnorm);
393:   VecDestroy(&r);

395:   PetscViewerPushFormat(viewer, format);
396:   PetscViewerASCIIAddTab(viewer, tablevel);
397:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
398:   PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " KSP %s resid norm %14.12e true resid norm %14.12e ||r(i)||/||b|| %14.12e\n", n, normtype, (double)rnorm, (double)truenorm, (double)(truenorm / bnorm));
399:   PetscViewerASCIISubtractTab(viewer, tablevel);
400:   PetscViewerPopFormat(viewer);
401:   return 0;
402: }

404: /*@C
405:   KSPMonitorTrueResidualDraw - Plots the true residual at each iteration of an iterative solver.

407:   Collective

409:   Input Parameters:
410: + ksp   - iterative context
411: . n     - iteration number
412: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
413: - vf    - The viewer context of type `PETSCVIEWERDRAW`

415:   Options Database Key:
416: . -ksp_monitor_true_residual draw - Activates `KSPMonitorResidualDraw()`

418:   Level: intermediate

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

424: .seealso: [](chapter_ksp), `PETSCVIEWERDRAW`, `KSP`, `KSPMonitorSet()`, `KSPMonitorTrueResidual()`, `KSPMonitorResidual()`, `KSPMonitorTrueResidualDrawLG()`
425: @*/
426: PetscErrorCode KSPMonitorTrueResidualDraw(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
427: {
428:   PetscViewer       viewer = vf->viewer;
429:   PetscViewerFormat format = vf->format;
430:   Vec               r;

433:   PetscViewerPushFormat(viewer, format);
434:   KSPBuildResidual(ksp, NULL, NULL, &r);
435:   PetscObjectSetName((PetscObject)r, "Residual");
436:   PetscObjectCompose((PetscObject)r, "__Vec_bc_zero__", (PetscObject)ksp);
437:   VecView(r, viewer);
438:   PetscObjectCompose((PetscObject)r, "__Vec_bc_zero__", NULL);
439:   VecDestroy(&r);
440:   PetscViewerPopFormat(viewer);
441:   return 0;
442: }

444: /*@C
445:   KSPMonitorTrueResidualDrawLG - Plots the true residual norm at each iteration of an iterative solver.

447:   Collective

449:   Input Parameters:
450: + ksp   - iterative context
451: . n     - iteration number
452: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
453: - vf    - The viewer context

455:   Options Database Key:
456: . -ksp_monitor_true_residual draw::draw_lg - Activates KSPMonitorTrueResidualDrawLG()

458:   Level: intermediate

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

464:   Call `KSPMonitorTrueResidualDrawLGCreate()` to create the context needed for this monitor

466: .seealso: [](chapter_ksp), `PETSCVIEWERDRAW`, `KSP`, `KSPMonitorSet()`, `KSPMonitorTrueResidual()`, `KSPMonitorTrueResidualDraw()`, `KSPMonitorResidual`,
467:           `KSPMonitorTrueResidualDrawLGCreate()`
468: @*/
469: PetscErrorCode KSPMonitorTrueResidualDrawLG(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
470: {
471:   PetscViewer        viewer = vf->viewer;
472:   PetscViewerFormat  format = vf->format;
473:   PetscDrawLG        lg     = vf->lg;
474:   Vec                r;
475:   KSPConvergedReason reason;
476:   PetscReal          truenorm, x[2], y[2];

480:   KSPBuildResidual(ksp, NULL, NULL, &r);
481:   VecNorm(r, NORM_2, &truenorm);
482:   VecDestroy(&r);
483:   PetscViewerPushFormat(viewer, format);
484:   if (!n) PetscDrawLGReset(lg);
485:   x[0] = (PetscReal)n;
486:   if (rnorm > 0.0) y[0] = PetscLog10Real(rnorm);
487:   else y[0] = -15.0;
488:   x[1] = (PetscReal)n;
489:   if (truenorm > 0.0) y[1] = PetscLog10Real(truenorm);
490:   else y[1] = -15.0;
491:   PetscDrawLGAddPoint(lg, x, y);
492:   KSPGetConvergedReason(ksp, &reason);
493:   if (n <= 20 || !(n % 5) || reason) {
494:     PetscDrawLGDraw(lg);
495:     PetscDrawLGSave(lg);
496:   }
497:   PetscViewerPopFormat(viewer);
498:   return 0;
499: }

501: /*@C
502:   KSPMonitorTrueResidualDrawLGCreate - Creates the context for the true residual monitor `KSPMonitorTrueResidualDrawLG()`

504:   Collective

506:   Input Parameters:
507: + viewer - The `PetscViewer` of type `PETSCVIEWERDRAW`
508: . format - The viewer format
509: - ctx    - An optional user context

511:   Output Parameter:
512: . vf    - The viewer context

514:   Level: intermediate

516: .seealso: [](chapter_ksp), `PETSCVIEWERDRAW`, `KSP`, `KSPMonitorSet()`, `KSPMonitorTrueResidual()`
517: @*/
518: PetscErrorCode KSPMonitorTrueResidualDrawLGCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
519: {
520:   const char *names[] = {"preconditioned", "true"};

522:   PetscViewerAndFormatCreate(viewer, format, vf);
523:   (*vf)->data = ctx;
524:   KSPMonitorLGCreate(PetscObjectComm((PetscObject)viewer), NULL, NULL, "Log Residual Norm", 2, names, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &(*vf)->lg);
525:   return 0;
526: }

528: /*@C
529:   KSPMonitorTrueResidualMax - Prints the true residual max norm at each iteration of an iterative solver.

531:   Collective

533:   Input Parameters:
534: + ksp   - iterative context
535: . n     - iteration number
536: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
537: - vf    - The viewer context

539:   Options Database Key:
540: . -ksp_monitor_true_residual_max - Activates `KSPMonitorTrueResidualMax()`

542:   Level: intermediate

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

548: .seealso: [](chapter_ksp), `KSP`, `KSPMonitorSet()`, `KSPMonitorResidual()`, `KSPMonitorTrueResidualMaxNorm()`
549: @*/
550: PetscErrorCode KSPMonitorTrueResidualMax(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
551: {
552:   PetscViewer       viewer = vf->viewer;
553:   PetscViewerFormat format = vf->format;
554:   Vec               r;
555:   PetscReal         truenorm, bnorm;
556:   char              normtype[256];
557:   PetscInt          tablevel;
558:   const char       *prefix;

561:   PetscObjectGetTabLevel((PetscObject)ksp, &tablevel);
562:   PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix);
563:   PetscStrncpy(normtype, KSPNormTypes[ksp->normtype], sizeof(normtype));
564:   PetscStrtolower(normtype);
565:   KSPBuildResidual(ksp, NULL, NULL, &r);
566:   VecNorm(r, NORM_INFINITY, &truenorm);
567:   VecNorm(ksp->vec_rhs, NORM_INFINITY, &bnorm);
568:   VecDestroy(&r);

570:   PetscViewerPushFormat(viewer, format);
571:   PetscViewerASCIIAddTab(viewer, tablevel);
572:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
573:   PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " KSP %s true resid max norm %14.12e ||r(i)||/||b|| %14.12e\n", n, normtype, (double)truenorm, (double)(truenorm / bnorm));
574:   PetscViewerASCIISubtractTab(viewer, tablevel);
575:   PetscViewerPopFormat(viewer);
576:   return 0;
577: }

579: /*@C
580:   KSPMonitorError - Prints the error norm, as well as the (possibly preconditioned) residual norm, at each iteration of an iterative solver.

582:   Collective

584:   Input Parameters:
585: + ksp   - iterative context
586: . n     - iteration number
587: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
588: - vf    - The viewer context

590:   Options Database Key:
591: . -ksp_monitor_error - Activates `KSPMonitorError()`

593:   Level: intermediate

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

599: .seealso: [](chapter_ksp), `KSP`, `KSPMonitorSet()`, `KSPMonitorResidual()`, `KSPMonitorTrueResidualMaxNorm()`
600: @*/
601: PetscErrorCode KSPMonitorError(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
602: {
603:   PetscViewer       viewer = vf->viewer;
604:   PetscViewerFormat format = vf->format;
605:   DM                dm;
606:   Vec               sol;
607:   PetscReal        *errors;
608:   PetscInt          Nf, f;
609:   PetscInt          tablevel;
610:   const char       *prefix;

613:   PetscObjectGetTabLevel((PetscObject)ksp, &tablevel);
614:   PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix);
615:   KSPGetDM(ksp, &dm);
616:   DMGetNumFields(dm, &Nf);
617:   DMGetGlobalVector(dm, &sol);
618:   KSPBuildSolution(ksp, sol, NULL);
619:   /* TODO: Make a different monitor that flips sign for SNES, Newton system is A dx = -b, so we need to negate the solution */
620:   VecScale(sol, -1.0);
621:   PetscCalloc1(Nf, &errors);
622:   DMComputeError(dm, sol, errors, NULL);

624:   PetscViewerPushFormat(viewer, format);
625:   PetscViewerASCIIAddTab(viewer, tablevel);
626:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Error norms for %s solve.\n", prefix);
627:   PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " KSP Error norm %s", n, Nf > 1 ? "[" : "");
628:   PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
629:   for (f = 0; f < Nf; ++f) {
630:     if (f > 0) PetscViewerASCIIPrintf(viewer, ", ");
631:     PetscViewerASCIIPrintf(viewer, "%14.12e", (double)errors[f]);
632:   }
633:   PetscViewerASCIIPrintf(viewer, "%s resid norm %14.12e\n", Nf > 1 ? "]" : "", (double)rnorm);
634:   PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
635:   PetscViewerASCIISubtractTab(viewer, tablevel);
636:   PetscViewerPopFormat(viewer);
637:   DMRestoreGlobalVector(dm, &sol);
638:   return 0;
639: }

641: /*@C
642:   KSPMonitorErrorDraw - Plots the error at each iteration of an iterative solver.

644:   Collective

646:   Input Parameters:
647: + ksp   - iterative context
648: . n     - iteration number
649: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
650: - vf    - The viewer context

652:   Options Database Key:
653: . -ksp_monitor_error draw - Activates `KSPMonitorErrorDraw()`

655:   Level: intermediate

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

661: .seealso: [](chapter_ksp), `PETSCVIEWERDRAW`, `KSP`, `KSPMonitorSet()`, `KSPMonitorTrueResidual()`, `KSPMonitorErrorDrawLG()`
662: @*/
663: PetscErrorCode KSPMonitorErrorDraw(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
664: {
665:   PetscViewer       viewer = vf->viewer;
666:   PetscViewerFormat format = vf->format;
667:   DM                dm;
668:   Vec               sol, e;

671:   PetscViewerPushFormat(viewer, format);
672:   KSPGetDM(ksp, &dm);
673:   DMGetGlobalVector(dm, &sol);
674:   KSPBuildSolution(ksp, sol, NULL);
675:   DMComputeError(dm, sol, NULL, &e);
676:   PetscObjectSetName((PetscObject)e, "Error");
677:   PetscObjectCompose((PetscObject)e, "__Vec_bc_zero__", (PetscObject)ksp);
678:   VecView(e, viewer);
679:   PetscObjectCompose((PetscObject)e, "__Vec_bc_zero__", NULL);
680:   VecDestroy(&e);
681:   DMRestoreGlobalVector(dm, &sol);
682:   PetscViewerPopFormat(viewer);
683:   return 0;
684: }

686: /*@C
687:   KSPMonitorErrorDrawLG - Plots the error and residual norm at each iteration of an iterative solver.

689:   Collective

691:   Input Parameters:
692: + ksp   - iterative context
693: . n     - iteration number
694: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
695: - vf    - The viewer context

697:   Options Database Key:
698: . -ksp_monitor_error draw::draw_lg - Activates `KSPMonitorTrueResidualDrawLG()`

700:   Level: intermediate

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

706:   Call `KSPMonitorErrorDrawLGCreate()` to create the context used with this monitor

708: .seealso: [](chapter_ksp), `PETSCVIEWERDRAW`, `KSP`, `KSPMonitorSet()`, `KSPMonitorTrueResidual()`, `KSPMonitorErrorDraw()`
709: @*/
710: PetscErrorCode KSPMonitorErrorDrawLG(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
711: {
712:   PetscViewer        viewer = vf->viewer;
713:   PetscViewerFormat  format = vf->format;
714:   PetscDrawLG        lg     = vf->lg;
715:   DM                 dm;
716:   Vec                sol;
717:   KSPConvergedReason reason;
718:   PetscReal         *x, *errors;
719:   PetscInt           Nf, f;

723:   KSPGetDM(ksp, &dm);
724:   DMGetNumFields(dm, &Nf);
725:   DMGetGlobalVector(dm, &sol);
726:   KSPBuildSolution(ksp, sol, NULL);
727:   /* TODO: Make a different monitor that flips sign for SNES, Newton system is A dx = -b, so we need to negate the solution */
728:   VecScale(sol, -1.0);
729:   PetscCalloc2(Nf + 1, &x, Nf + 1, &errors);
730:   DMComputeError(dm, sol, errors, NULL);

732:   PetscViewerPushFormat(viewer, format);
733:   if (!n) PetscDrawLGReset(lg);
734:   for (f = 0; f < Nf; ++f) {
735:     x[f]      = (PetscReal)n;
736:     errors[f] = errors[f] > 0.0 ? PetscLog10Real(errors[f]) : -15.;
737:   }
738:   x[Nf]      = (PetscReal)n;
739:   errors[Nf] = rnorm > 0.0 ? PetscLog10Real(rnorm) : -15.;
740:   PetscDrawLGAddPoint(lg, x, errors);
741:   KSPGetConvergedReason(ksp, &reason);
742:   if (n <= 20 || !(n % 5) || reason) {
743:     PetscDrawLGDraw(lg);
744:     PetscDrawLGSave(lg);
745:   }
746:   PetscViewerPopFormat(viewer);
747:   return 0;
748: }

750: /*@C
751:   KSPMonitorErrorDrawLGCreate - Creates the context for the error and preconditioned residual plotter `KSPMonitorErrorDrawLG()`

753:   Collective

755:   Input Parameters:
756: + viewer - The `PetscViewer`
757: . format - The viewer format
758: - ctx    - An optional user context

760:   Output Parameter:
761: . vf    - The viewer context

763:   Level: intermediate

765: .seealso: [](chapter_ksp), `PETSCVIEWERDRAW`, `KSP`, `KSPMonitorSet()`, `KSPMonitorTrueResidual()`, `KSPMonitorErrorDrawLG()`, `KSPMonitorErrorDrawLG()`
766: @*/
767: PetscErrorCode KSPMonitorErrorDrawLGCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
768: {
769:   KSP      ksp = (KSP)ctx;
770:   DM       dm;
771:   char   **names;
772:   PetscInt Nf, f;

774:   KSPGetDM(ksp, &dm);
775:   DMGetNumFields(dm, &Nf);
776:   PetscMalloc1(Nf + 1, &names);
777:   for (f = 0; f < Nf; ++f) {
778:     PetscObject disc;
779:     const char *fname;
780:     char        lname[PETSC_MAX_PATH_LEN];

782:     DMGetField(dm, f, NULL, &disc);
783:     PetscObjectGetName(disc, &fname);
784:     PetscStrncpy(lname, fname, PETSC_MAX_PATH_LEN);
785:     PetscStrlcat(lname, " Error", PETSC_MAX_PATH_LEN);
786:     PetscStrallocpy(lname, &names[f]);
787:   }
788:   PetscStrallocpy("residual", &names[Nf]);
789:   PetscViewerAndFormatCreate(viewer, format, vf);
790:   (*vf)->data = ctx;
791:   KSPMonitorLGCreate(PetscObjectComm((PetscObject)viewer), NULL, NULL, "Log Error Norm", Nf + 1, (const char **)names, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &(*vf)->lg);
792:   for (f = 0; f <= Nf; ++f) PetscFree(names[f]);
793:   PetscFree(names);
794:   return 0;
795: }

797: /*@C
798:   KSPMonitorSolution - Print the solution norm at each iteration of an iterative solver.

800:   Collective

802:   Input Parameters:
803: + ksp   - iterative context
804: . n     - iteration number
805: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
806: - vf    - The viewer context

808:   Options Database Key:
809: . -ksp_monitor_solution - Activates `KSPMonitorSolution()`

811:   Level: intermediate

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

817: .seealso: [](chapter_ksp), `KSPMonitorSet()`, `KSPMonitorTrueResidual()`
818: @*/
819: PetscErrorCode KSPMonitorSolution(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
820: {
821:   PetscViewer       viewer = vf->viewer;
822:   PetscViewerFormat format = vf->format;
823:   Vec               x;
824:   PetscReal         snorm;
825:   PetscInt          tablevel;
826:   const char       *prefix;

829:   KSPBuildSolution(ksp, NULL, &x);
830:   VecNorm(x, NORM_2, &snorm);
831:   PetscObjectGetTabLevel((PetscObject)ksp, &tablevel);
832:   PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix);
833:   PetscViewerPushFormat(viewer, format);
834:   PetscViewerASCIIAddTab(viewer, tablevel);
835:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Solution norms for %s solve.\n", prefix);
836:   PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " KSP Solution norm %14.12e \n", n, (double)snorm);
837:   PetscViewerASCIISubtractTab(viewer, tablevel);
838:   PetscViewerPopFormat(viewer);
839:   return 0;
840: }

842: /*@C
843:   KSPMonitorSolutionDraw - Plots the solution at each iteration of an iterative solver.

845:   Collective

847:   Input Parameters:
848: + ksp   - iterative context
849: . n     - iteration number
850: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
851: - vf    - The viewer context

853:   Options Database Key:
854: . -ksp_monitor_solution draw - Activates `KSPMonitorSolutionDraw()`

856:   Level: intermediate

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

862: .seealso: [](chapter_ksp), `KSPMonitorSet()`, `KSPMonitorTrueResidual()`
863: @*/
864: PetscErrorCode KSPMonitorSolutionDraw(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
865: {
866:   PetscViewer       viewer = vf->viewer;
867:   PetscViewerFormat format = vf->format;
868:   Vec               x;

871:   KSPBuildSolution(ksp, NULL, &x);
872:   PetscViewerPushFormat(viewer, format);
873:   PetscObjectSetName((PetscObject)x, "Solution");
874:   PetscObjectCompose((PetscObject)x, "__Vec_bc_zero__", (PetscObject)ksp);
875:   VecView(x, viewer);
876:   PetscObjectCompose((PetscObject)x, "__Vec_bc_zero__", NULL);
877:   PetscViewerPopFormat(viewer);
878:   return 0;
879: }

881: /*@C
882:   KSPMonitorSolutionDrawLG - Plots the solution norm at each iteration of an iterative solver.

884:   Collective

886:   Input Parameters:
887: + ksp   - iterative context
888: . n     - iteration number
889: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
890: - vf    - The viewer context

892:   Options Database Key:
893: . -ksp_monitor_solution draw::draw_lg - Activates `KSPMonitorSolutionDrawLG()`

895:   Level: intermediate

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

901:    Call `KSPMonitorSolutionDrawLGCreate()` to create the context needed with this monitor

903: .seealso: [](chapter_ksp), `KSP`, `KSPMonitorSet()`, `KSPMonitorTrueResidual()`, `KSPMonitorSolutionDrawLGCreate()`
904: @*/
905: PetscErrorCode KSPMonitorSolutionDrawLG(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
906: {
907:   PetscViewer        viewer = vf->viewer;
908:   PetscViewerFormat  format = vf->format;
909:   PetscDrawLG        lg     = vf->lg;
910:   Vec                u;
911:   KSPConvergedReason reason;
912:   PetscReal          snorm, x, y;

916:   KSPBuildSolution(ksp, NULL, &u);
917:   VecNorm(u, NORM_2, &snorm);
918:   PetscViewerPushFormat(viewer, format);
919:   if (!n) PetscDrawLGReset(lg);
920:   x = (PetscReal)n;
921:   if (snorm > 0.0) y = PetscLog10Real(snorm);
922:   else y = -15.0;
923:   PetscDrawLGAddPoint(lg, &x, &y);
924:   KSPGetConvergedReason(ksp, &reason);
925:   if (n <= 20 || !(n % 5) || reason) {
926:     PetscDrawLGDraw(lg);
927:     PetscDrawLGSave(lg);
928:   }
929:   PetscViewerPopFormat(viewer);
930:   return 0;
931: }

933: /*@C
934:   KSPMonitorSolutionDrawLGCreate - Creates the context for the `KSP` monitor `KSPMonitorSolutionDrawLG()`

936:   Collective

938:   Input Parameters:
939: + viewer - The `PetscViewer`
940: . format - The viewer format
941: - ctx    - An optional user context

943:   Output Parameter:
944: . vf    - The viewer context

946:   Level: intermediate

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

952: .seealso: [](chapter_ksp), `KSPMonitorSet()`, `KSPMonitorTrueResidual()`
953: @*/
954: PetscErrorCode KSPMonitorSolutionDrawLGCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
955: {
956:   PetscViewerAndFormatCreate(viewer, format, vf);
957:   (*vf)->data = ctx;
958:   KSPMonitorLGCreate(PetscObjectComm((PetscObject)viewer), NULL, NULL, "Log Solution Norm", 1, NULL, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &(*vf)->lg);
959:   return 0;
960: }

962: /*@C
963:   KSPMonitorSingularValue - Prints the two norm of the true residual and estimation of the extreme singular values of the preconditioned problem at each iteration.

965:   Logically Collective

967:   Input Parameters:
968: + ksp   - the iterative context
969: . n     - the iteration
970: . rnorm - the two norm of the residual
971: - vf    - The viewer context

973:   Options Database Key:
974: . -ksp_monitor_singular_value - Activates `KSPMonitorSingularValue()`

976:   Level: intermediate

978:   Notes:
979:   The `KSPCG` solver uses the Lanczos technique for eigenvalue computation,
980:   while `KSPGMRES` uses the Arnoldi technique; other iterative methods do
981:   not currently compute singular values.

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

986:   Call `KSPMonitorSingularValueCreate()` to create the context needed by this monitor

988: .seealso: [](chapter_ksp), `KSP`, `KSPMonitorSet()`, `KSPComputeExtremeSingularValues()`, `KSPMonitorSingularValueCreate()`
989: @*/
990: PetscErrorCode KSPMonitorSingularValue(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
991: {
992:   PetscViewer       viewer = vf->viewer;
993:   PetscViewerFormat format = vf->format;
994:   PetscReal         emin, emax;
995:   PetscInt          tablevel;
996:   const char       *prefix;

1000:   PetscObjectGetTabLevel((PetscObject)ksp, &tablevel);
1001:   PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix);
1002:   PetscViewerPushFormat(viewer, format);
1003:   PetscViewerASCIIAddTab(viewer, tablevel);
1004:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
1005:   if (!ksp->calc_sings) {
1006:     PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " KSP Residual norm %14.12e \n", n, (double)rnorm);
1007:   } else {
1008:     KSPComputeExtremeSingularValues(ksp, &emax, &emin);
1009:     PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " KSP Residual norm %14.12e %% max %14.12e min %14.12e max/min %14.12e\n", n, (double)rnorm, (double)emax, (double)emin, (double)(emax / emin));
1010:   }
1011:   PetscViewerASCIISubtractTab(viewer, tablevel);
1012:   PetscViewerPopFormat(viewer);
1013:   return 0;
1014: }

1016: /*@C
1017:   KSPMonitorSingularValueCreate - Creates the singular value monitor context needed by `KSPMonitorSingularValue()`

1019:   Collective

1021:   Input Parameters:
1022: + viewer - The PetscViewer
1023: . format - The viewer format
1024: - ctx    - An optional user context

1026:   Output Parameter:
1027: . vf    - The viewer context

1029:   Level: intermediate

1031: .seealso: [](chapter_ksp), `KSP`, `KSPMonitorSet()`, `KSPMonitorSingularValue()`
1032: @*/
1033: PetscErrorCode KSPMonitorSingularValueCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
1034: {
1035:   KSP ksp = (KSP)ctx;

1037:   PetscViewerAndFormatCreate(viewer, format, vf);
1038:   (*vf)->data = ctx;
1039:   KSPSetComputeSingularValues(ksp, PETSC_TRUE);
1040:   return 0;
1041: }

1043: /*@C
1044:    KSPMonitorDynamicTolerance - Recompute the inner tolerance in every outer iteration in an adaptive way.

1046:    Collective

1048:    Input Parameters:
1049: +  ksp   - iterative context
1050: .  n     - iteration number (not used)
1051: .  fnorm - the current residual norm
1052: -  dummy - some context as a C struct. fields:
1053:              coef: a scaling coefficient. default 1.0. can be passed through
1054:                    -sub_ksp_dynamic_tolerance_param
1055:              bnrm: norm of the right-hand side. store it to avoid repeated calculation

1057:    Notes:
1058:    This may be useful for a flexibly preconditioner Krylov method to
1059:    control the accuracy of the inner solves needed to guarantee the
1060:    convergence of the outer iterations.

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

1065:    Level: advanced

1067: .seealso: `KSPMonitorDynamicToleranceDestroy()`
1068: @*/
1069: PetscErrorCode KSPMonitorDynamicTolerance(KSP ksp, PetscInt its, PetscReal fnorm, void *dummy)
1070: {
1071:   PC            pc;
1072:   PetscReal     outer_rtol, outer_abstol, outer_dtol, inner_rtol;
1073:   PetscInt      outer_maxits, nksp, first, i;
1074:   KSPDynTolCtx *scale  = (KSPDynTolCtx *)dummy;
1075:   KSP          *subksp = NULL;
1076:   KSP           kspinner;
1077:   PetscBool     flg;

1079:   KSPGetPC(ksp, &pc);

1081:   /* compute inner_rtol */
1082:   if (scale->bnrm < 0.0) {
1083:     Vec b;
1084:     KSPGetRhs(ksp, &b);
1085:     VecNorm(b, NORM_2, &(scale->bnrm));
1086:   }
1087:   KSPGetTolerances(ksp, &outer_rtol, &outer_abstol, &outer_dtol, &outer_maxits);
1088:   inner_rtol = PetscMin(scale->coef * scale->bnrm * outer_rtol / fnorm, 0.999);
1089:   /* PetscCall(PetscPrintf(PETSC_COMM_WORLD, "        Inner rtol = %g\n",
1090:      (double)inner_rtol)); */

1092:   /* if pc is ksp */
1093:   PetscObjectTypeCompare((PetscObject)pc, PCKSP, &flg);
1094:   if (flg) {
1095:     PCKSPGetKSP(pc, &kspinner);
1096:     KSPSetTolerances(kspinner, inner_rtol, outer_abstol, outer_dtol, outer_maxits);
1097:     return 0;
1098:   }

1100:   /* if pc is bjacobi */
1101:   PetscObjectTypeCompare((PetscObject)pc, PCBJACOBI, &flg);
1102:   if (flg) {
1103:     PCBJacobiGetSubKSP(pc, &nksp, &first, &subksp);
1104:     if (subksp) {
1105:       for (i = 0; i < nksp; i++) KSPSetTolerances(subksp[i], inner_rtol, outer_abstol, outer_dtol, outer_maxits);
1106:       return 0;
1107:     }
1108:   }

1110:   /* if pc is deflation*/
1111:   PetscObjectTypeCompare((PetscObject)pc, PCDEFLATION, &flg);
1112:   if (flg) {
1113:     PCDeflationGetCoarseKSP(pc, &kspinner);
1114:     KSPSetTolerances(kspinner, inner_rtol, outer_abstol, outer_dtol, PETSC_DEFAULT);
1115:     return 0;
1116:   }

1118:   /* todo: dynamic tolerance may apply to other types of pc too */
1119:   return 0;
1120: }

1122: /*
1123:   Destroy the dummy context used in KSPMonitorDynamicTolerance()
1124: */
1125: PetscErrorCode KSPMonitorDynamicToleranceDestroy(void **dummy)
1126: {
1127:   PetscFree(*dummy);
1128:   return 0;
1129: }

1131: /*@C
1132:    KSPConvergedSkip - Convergence test that do not return as converged
1133:    until the maximum number of iterations is reached.

1135:    Collective

1137:    Input Parameters:
1138: +  ksp   - iterative context
1139: .  n     - iteration number
1140: .  rnorm - 2-norm residual value (may be estimated)
1141: -  dummy - unused convergence context

1143:    Returns:
1144: .  reason - `KSP_CONVERGED_ITERATING`, `KSP_CONVERGED_ITS`

1146:    Notes:
1147:    This should be used as the convergence test with the option
1148:    `KSPSetNormType`(ksp,`KSP_NORM_NONE`), since norms of the residual are
1149:    not computed. Convergence is then declared after the maximum number
1150:    of iterations have been reached. Useful when one is using `KSPCG` or
1151:    `KSPBCGS`. [](sec_flexibleksp)

1153:    Level: advanced

1155: .seealso: [](chapter_ksp), `KSP`, `KSPCG`, `KSPBCGS`, `KSPSetConvergenceTest()`, `KSPSetTolerances()`, `KSPSetNormType()`, [](sec_flexibleksp),
1156: @*/
1157: PetscErrorCode KSPConvergedSkip(KSP ksp, PetscInt n, PetscReal rnorm, KSPConvergedReason *reason, void *dummy)
1158: {
1161:   *reason = KSP_CONVERGED_ITERATING;
1162:   if (n >= ksp->max_it) *reason = KSP_CONVERGED_ITS;
1163:   return 0;
1164: }

1166: /*@C
1167:    KSPConvergedDefaultCreate - Creates and initializes the context used by the `KSPConvergedDefault()` function

1169:    Note Collective

1171:    Output Parameter:
1172: .  ctx - convergence context

1174:    Level: intermediate

1176: .seealso: [](chapter_ksp), `KSP`, `KSPConvergedDefault()`, `KSPConvergedDefaultDestroy()`, `KSPSetConvergenceTest()`, `KSPSetTolerances()`,
1177:           `KSPConvergedSkip()`, `KSPConvergedReason`, `KSPGetConvergedReason()`, `KSPConvergedDefaultSetUIRNorm()`, `KSPConvergedDefaultSetUMIRNorm()`, `KSPConvergedDefaultSetConvergedMaxits()`
1178: @*/
1179: PetscErrorCode KSPConvergedDefaultCreate(void **ctx)
1180: {
1181:   KSPConvergedDefaultCtx *cctx;

1183:   PetscNew(&cctx);
1184:   *ctx = cctx;
1185:   return 0;
1186: }

1188: /*@
1189:    KSPConvergedDefaultSetUIRNorm - makes the default convergence test use || B*(b - A*(initial guess))||
1190:       instead of || B*b ||. In the case of right preconditioner or if `KSPSetNormType`(ksp,`KSP_NORM_UNPRECONDITIONED`)
1191:       is used there is no B in the above formula. UIRNorm is short for Use Initial Residual Norm.

1193:    Collective

1195:    Input Parameters:
1196: .  ksp   - iterative context

1198:    Options Database Key:
1199: .   -ksp_converged_use_initial_residual_norm <bool> - Use initial residual norm for computing relative convergence

1201:    Notes:
1202:    Use `KSPSetTolerances()` to alter the defaults for rtol, abstol, dtol.

1204:    The precise values of reason are macros such as `KSP_CONVERGED_RTOL`, which
1205:    are defined in petscksp.h.

1207:    If the convergence test is not `KSPConvergedDefault()` then this is ignored.

1209:    If right preconditioning is being used then B does not appear in the above formula.

1211:    Level: intermediate

1213: .seealso: [](chapter_ksp), `KSP`, `KSPSetConvergenceTest()`, `KSPSetTolerances()`, `KSPConvergedSkip()`, `KSPConvergedReason`, `KSPGetConvergedReason()`, `KSPConvergedDefaultSetUMIRNorm()`, `KSPConvergedDefaultSetConvergedMaxits()`
1214: @*/
1215: PetscErrorCode KSPConvergedDefaultSetUIRNorm(KSP ksp)
1216: {
1217:   KSPConvergedDefaultCtx *ctx = (KSPConvergedDefaultCtx *)ksp->cnvP;

1220:   if (ksp->converged != KSPConvergedDefault) return 0;
1222:   ctx->initialrtol = PETSC_TRUE;
1223:   return 0;
1224: }

1226: /*@
1227:    KSPConvergedDefaultSetUMIRNorm - makes the default convergence test use min(|| B*(b - A*(initial guess))||,|| B*b ||)
1228:       In the case of right preconditioner or if `KSPSetNormType`(ksp,`KSP_NORM_UNPRECONDITIONED`)
1229:       is used there is no B in the above formula. UMIRNorm is short for Use Minimum Initial Residual Norm.

1231:    Collective

1233:    Input Parameters:
1234: .  ksp   - iterative context

1236:    Options Database Key:
1237: .   -ksp_converged_use_min_initial_residual_norm <bool> - Use minimum of initial residual norm and b for computing relative convergence

1239:    Use `KSPSetTolerances()` to alter the defaults for rtol, abstol, dtol.

1241:    The precise values of reason are macros such as `KSP_CONVERGED_RTOL`, which
1242:    are defined in petscksp.h.

1244:    Level: intermediate

1246: .seealso: [](chapter_ksp), `KSP`, `KSPSetConvergenceTest()`, `KSPSetTolerances()`, `KSPConvergedSkip()`, `KSPConvergedReason`, `KSPGetConvergedReason()`, `KSPConvergedDefaultSetUIRNorm()`, `KSPConvergedDefaultSetConvergedMaxits()`
1247: @*/
1248: PetscErrorCode KSPConvergedDefaultSetUMIRNorm(KSP ksp)
1249: {
1250:   KSPConvergedDefaultCtx *ctx = (KSPConvergedDefaultCtx *)ksp->cnvP;

1253:   if (ksp->converged != KSPConvergedDefault) return 0;
1255:   ctx->mininitialrtol = PETSC_TRUE;
1256:   return 0;
1257: }

1259: /*@
1260:    KSPConvergedDefaultSetConvergedMaxits - allows the default convergence test to declare convergence and return `KSP_CONVERGED_ITS` if the maximum number of iterations is reached

1262:    Collective

1264:    Input Parameters:
1265: +  ksp - iterative context
1266: -  flg - boolean flag

1268:    Options Database Key:
1269: .   -ksp_converged_maxits <bool> - Declare convergence if the maximum number of iterations is reached

1271:    The precise values of reason are macros such as `KSP_CONVERGED_RTOL`, which
1272:    are defined in petscksp.h.

1274:    Level: intermediate

1276: .seealso: [](chapter_ksp), `KSP`, `KSPSetConvergenceTest()`, `KSPSetTolerances()`, `KSPConvergedSkip()`, `KSPConvergedReason`, `KSPGetConvergedReason()`, `KSPConvergedDefaultSetUMIRNorm()`, `KSPConvergedDefaultSetUIRNorm()`
1277: @*/
1278: PetscErrorCode KSPConvergedDefaultSetConvergedMaxits(KSP ksp, PetscBool flg)
1279: {
1280:   KSPConvergedDefaultCtx *ctx = (KSPConvergedDefaultCtx *)ksp->cnvP;

1284:   if (ksp->converged != KSPConvergedDefault) return 0;
1285:   ctx->convmaxits = flg;
1286:   return 0;
1287: }

1289: /*@C
1290:    KSPConvergedDefault - Default code to determine convergence of the linear iterative solvers

1292:    Collective

1294:    Input Parameters:
1295: +  ksp   - iterative context
1296: .  n     - iteration number
1297: .  rnorm - residual norm (may be estimated, depending on the method may be the preconditioned residual norm)
1298: -  ctx - convergence context which must be created by `KSPConvergedDefaultCreate()`

1300:    Output Parameter:
1301: .  reason - the convergence reason; it is positive if the iteration has converged,
1302:             negative if the iteration has diverged, and `KSP_CONVERGED_ITERATING` otherwise

1304:    Options Database Keys:
1305: +   -ksp_max_it - maximum number of linear iterations
1306: .   -ksp_rtol rtol - relative tolerance used in default determination of convergence, i.e. if residual norm decreases by this factor than convergence is declared
1307: .   -ksp_atol abstol - absolute tolerance used in default convergence test, i.e. if residual norm is less than this then convergence is declared
1308: .   -ksp_divtol tol - if residual norm increases by this factor than divergence is declared
1309: .   -ksp_converged_use_initial_residual_norm - see `KSPConvergedDefaultSetUIRNorm()`
1310: .   -ksp_converged_use_min_initial_residual_norm - see `KSPConvergedDefaultSetUMIRNorm()`
1311: -   -ksp_converged_maxits - see `KSPConvergedDefaultSetConvergedMaxits()`

1313:    Level: advanced

1315:    Notes:
1316:    `KSPConvergedDefault()` reaches convergence when   rnorm < MAX (rtol * rnorm_0, abstol);
1317:    Divergence is detected if rnorm > dtol * rnorm_0, or when failures are detected throughout the iteration.
1318:    By default, reaching the maximum number of iterations is considered divergence (i.e. KSP_DIVERGED_ITS).
1319:    In order to have PETSc declaring convergence in such a case (i.e. `KSP_CONVERGED_ITS`), users can use `KSPConvergedDefaultSetConvergedMaxits()`

1321:    where:
1322: +     rtol - relative tolerance,
1323: .     abstol - absolute tolerance.
1324: .     dtol - divergence tolerance,
1325: -     rnorm_0 - the two norm of the right hand side (or the preconditioned norm, depending on what was set with
1326:           `KSPSetNormType()`. When initial guess is non-zero you
1327:           can call `KSPConvergedDefaultSetUIRNorm()` to use the norm of (b - A*(initial guess))
1328:           as the starting point for relative norm convergence testing, that is as rnorm_0.
1329:           Call `KSPConvergedDefaultSetUMIRNorm()` to use the minimum of the norm of (b - A*(initial guess)) and the norm of b as the starting point.

1331:    Use `KSPSetTolerances()` to alter the defaults for rtol, abstol, dtol.

1333:    Use `KSPSetNormType()` (or -ksp_norm_type <none,preconditioned,unpreconditioned,natural>) to change the norm used for computing rnorm

1335:    The precise values of reason are macros such as `KSP_CONVERGED_RTOL`, which are defined in petscksp.h.

1337:    This routine is used by `KSP` by default so the user generally never needs call it directly.

1339:    Use `KSPSetConvergenceTest()` to provide your own test instead of using this one.

1341:    Call `KSPSetConvergenceTest()` with the ctx, as created above and the destruction function `KSPConvergedDefaultDestroy()`

1343: .seealso: [](chapter_ksp), `KSP`, `KSPSetConvergenceTest()`, `KSPSetTolerances()`, `KSPConvergedSkip()`, `KSPConvergedReason`, `KSPGetConvergedReason()`,
1344:           `KSPConvergedDefaultSetUIRNorm()`, `KSPConvergedDefaultSetUMIRNorm()`, `KSPConvergedDefaultSetConvergedMaxits()`, `KSPConvergedDefaultCreate()`, `KSPConvergedDefaultDestroy()`
1345: @*/
1346: PetscErrorCode KSPConvergedDefault(KSP ksp, PetscInt n, PetscReal rnorm, KSPConvergedReason *reason, void *ctx)
1347: {
1348:   KSPConvergedDefaultCtx *cctx = (KSPConvergedDefaultCtx *)ctx;
1349:   KSPNormType             normtype;

1355:   *reason = KSP_CONVERGED_ITERATING;

1357:   if (cctx->convmaxits && n >= ksp->max_it) {
1358:     *reason = KSP_CONVERGED_ITS;
1359:     PetscInfo(ksp, "Linear solver has converged. Maximum number of iterations reached %" PetscInt_FMT "\n", n);
1360:     return 0;
1361:   }
1362:   KSPGetNormType(ksp, &normtype);
1363:   if (normtype == KSP_NORM_NONE) return 0;

1365:   if (!n) {
1366:     /* if user gives initial guess need to compute norm of b */
1367:     if (!ksp->guess_zero && !cctx->initialrtol) {
1368:       PetscReal snorm = 0.0;
1369:       if (ksp->normtype == KSP_NORM_UNPRECONDITIONED || ksp->pc_side == PC_RIGHT) {
1370:         PetscInfo(ksp, "user has provided nonzero initial guess, computing 2-norm of RHS\n");
1371:         VecNorm(ksp->vec_rhs, NORM_2, &snorm); /*     <- b'*b */
1372:       } else {
1373:         Vec z;
1374:         /* Should avoid allocating the z vector each time but cannot stash it in cctx because if KSPReset() is called the vector size might change */
1375:         VecDuplicate(ksp->vec_rhs, &z);
1376:         KSP_PCApply(ksp, ksp->vec_rhs, z);
1377:         if (ksp->normtype == KSP_NORM_PRECONDITIONED) {
1378:           PetscInfo(ksp, "user has provided nonzero initial guess, computing 2-norm of preconditioned RHS\n");
1379:           VecNorm(z, NORM_2, &snorm); /*    dp <- b'*B'*B*b */
1380:         } else if (ksp->normtype == KSP_NORM_NATURAL) {
1381:           PetscScalar norm;
1382:           PetscInfo(ksp, "user has provided nonzero initial guess, computing natural norm of RHS\n");
1383:           VecDot(ksp->vec_rhs, z, &norm);
1384:           snorm = PetscSqrtReal(PetscAbsScalar(norm)); /*    dp <- b'*B*b */
1385:         }
1386:         VecDestroy(&z);
1387:       }
1388:       /* handle special case of zero RHS and nonzero guess */
1389:       if (!snorm) {
1390:         PetscInfo(ksp, "Special case, user has provided nonzero initial guess and zero RHS\n");
1391:         snorm = rnorm;
1392:       }
1393:       if (cctx->mininitialrtol) ksp->rnorm0 = PetscMin(snorm, rnorm);
1394:       else ksp->rnorm0 = snorm;
1395:     } else {
1396:       ksp->rnorm0 = rnorm;
1397:     }
1398:     ksp->ttol = PetscMax(ksp->rtol * ksp->rnorm0, ksp->abstol);
1399:   }

1401:   if (n <= ksp->chknorm) return 0;

1403:   if (PetscIsInfOrNanReal(rnorm)) {
1404:     PCFailedReason pcreason;
1405:     PetscInt       sendbuf, recvbuf;
1406:     PCGetFailedReasonRank(ksp->pc, &pcreason);
1407:     sendbuf = (PetscInt)pcreason;
1408:     MPI_Allreduce(&sendbuf, &recvbuf, 1, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject)ksp));
1409:     if (recvbuf) {
1410:       *reason = KSP_DIVERGED_PC_FAILED;
1411:       PCSetFailedReason(ksp->pc, (PCFailedReason)recvbuf);
1412:       PetscInfo(ksp, "Linear solver pcsetup fails, declaring divergence \n");
1413:     } else {
1414:       *reason = KSP_DIVERGED_NANORINF;
1415:       PetscInfo(ksp, "Linear solver has created a not a number (NaN) as the residual norm, declaring divergence \n");
1416:     }
1417:   } else if (rnorm <= ksp->ttol) {
1418:     if (rnorm < ksp->abstol) {
1419:       PetscInfo(ksp, "Linear solver has converged. Residual norm %14.12e is less than absolute tolerance %14.12e at iteration %" PetscInt_FMT "\n", (double)rnorm, (double)ksp->abstol, n);
1420:       *reason = KSP_CONVERGED_ATOL;
1421:     } else {
1422:       if (cctx->initialrtol) {
1423:         PetscInfo(ksp, "Linear solver has converged. Residual norm %14.12e is less than relative tolerance %14.12e times initial residual norm %14.12e at iteration %" PetscInt_FMT "\n", (double)rnorm, (double)ksp->rtol, (double)ksp->rnorm0, n);
1424:       } else {
1425:         PetscInfo(ksp, "Linear solver has converged. Residual norm %14.12e is less than relative tolerance %14.12e times initial right hand side norm %14.12e at iteration %" PetscInt_FMT "\n", (double)rnorm, (double)ksp->rtol, (double)ksp->rnorm0, n);
1426:       }
1427:       *reason = KSP_CONVERGED_RTOL;
1428:     }
1429:   } else if (rnorm >= ksp->divtol * ksp->rnorm0) {
1430:     PetscInfo(ksp, "Linear solver is diverging. Initial right hand size norm %14.12e, current residual norm %14.12e at iteration %" PetscInt_FMT "\n", (double)ksp->rnorm0, (double)rnorm, n);
1431:     *reason = KSP_DIVERGED_DTOL;
1432:   }
1433:   return 0;
1434: }

1436: /*@C
1437:    KSPConvergedDefaultDestroy - Frees the space used by the `KSPConvergedDefault()` function context

1439:    Not Collective

1441:    Input Parameters:
1442: .  ctx - convergence context

1444:    Level: intermediate

1446:    Note:
1447:    Pass this function name into `KSPSetConvergenceTest()` along with the context obtained with `KSPConvergedDefaultCreate()` and `KSPConvergedDefault()`

1449: .seealso: [](chapter_ksp), `KSP`, `KSPConvergedDefault()`, `KSPConvergedDefaultCreate()`, `KSPSetConvergenceTest()`, `KSPSetTolerances()`, `KSPConvergedSkip()`,
1450:           `KSPConvergedReason`, `KSPGetConvergedReason()`, `KSPConvergedDefaultSetUIRNorm()`, `KSPConvergedDefaultSetUMIRNorm()`
1451: @*/
1452: PetscErrorCode KSPConvergedDefaultDestroy(void *ctx)
1453: {
1454:   KSPConvergedDefaultCtx *cctx = (KSPConvergedDefaultCtx *)ctx;

1456:   VecDestroy(&cctx->work);
1457:   PetscFree(ctx);
1458:   return 0;
1459: }

1461: /*
1462:    KSPBuildSolutionDefault - Default code to create/move the solution.

1464:    Collective

1466:    Input Parameters:
1467: +  ksp - iterative context
1468: -  v   - pointer to the user's vector

1470:    Output Parameter:
1471: .  V - pointer to a vector containing the solution

1473:    Level: advanced

1475:    Developers Note:
1476:    This is PETSC_EXTERN because it may be used by user written plugin KSP implementations

1478: .seealso: [](chapter_ksp), `KSPGetSolution()`, `KSPBuildResidualDefault()`
1479: */
1480: PetscErrorCode KSPBuildSolutionDefault(KSP ksp, Vec v, Vec *V)
1481: {
1482:   if (ksp->pc_side == PC_RIGHT) {
1483:     if (ksp->pc) {
1484:       if (v) {
1485:         KSP_PCApply(ksp, ksp->vec_sol, v);
1486:         *V = v;
1487:       } else SETERRQ(PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "Not working with right preconditioner");
1488:     } else {
1489:       if (v) {
1490:         VecCopy(ksp->vec_sol, v);
1491:         *V = v;
1492:       } else *V = ksp->vec_sol;
1493:     }
1494:   } else if (ksp->pc_side == PC_SYMMETRIC) {
1495:     if (ksp->pc) {
1497:       if (v) {
1498:         PCApplySymmetricRight(ksp->pc, ksp->vec_sol, v);
1499:         *V = v;
1500:       } else SETERRQ(PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "Not working with symmetric preconditioner");
1501:     } else {
1502:       if (v) {
1503:         VecCopy(ksp->vec_sol, v);
1504:         *V = v;
1505:       } else *V = ksp->vec_sol;
1506:     }
1507:   } else {
1508:     if (v) {
1509:       VecCopy(ksp->vec_sol, v);
1510:       *V = v;
1511:     } else *V = ksp->vec_sol;
1512:   }
1513:   return 0;
1514: }

1516: /*
1517:    KSPBuildResidualDefault - Default code to compute the residual.

1519:    Collecive on ksp

1521:    Input Parameters:
1522: .  ksp - iterative context
1523: .  t   - pointer to temporary vector
1524: .  v   - pointer to user vector

1526:    Output Parameter:
1527: .  V - pointer to a vector containing the residual

1529:    Level: advanced

1531:    Developers Note:
1532:    This is PETSC_EXTERN because it may be used by user written plugin KSP implementations

1534: .seealso: [](chapter_ksp), `KSPBuildSolutionDefault()`
1535: */
1536: PetscErrorCode KSPBuildResidualDefault(KSP ksp, Vec t, Vec v, Vec *V)
1537: {
1538:   Mat Amat, Pmat;

1540:   if (!ksp->pc) KSPGetPC(ksp, &ksp->pc);
1541:   PCGetOperators(ksp->pc, &Amat, &Pmat);
1542:   KSPBuildSolution(ksp, t, NULL);
1543:   KSP_MatMult(ksp, Amat, t, v);
1544:   VecAYPX(v, -1.0, ksp->vec_rhs);
1545:   *V = v;
1546:   return 0;
1547: }

1549: /*@C
1550:   KSPCreateVecs - Gets a number of work vectors suitably sized for the operator in the `KSP`

1552:   Collective

1554:   Input Parameters:
1555: + ksp  - iterative context
1556: . rightn  - number of right work vectors
1557: - leftn   - number of left work vectors to allocate

1559:   Output Parameters:
1560: +  right - the array of vectors created
1561: -  left - the array of left vectors

1563:    Level: advanced

1565:    Notes:
1566:    The right vector has as many elements as the matrix has columns. The left
1567:      vector has as many elements as the matrix has rows.

1569:    The vectors are new vectors that are not owned by the `KSP`, they should be destroyed with calls to `VecDestroyVecs()` when no longer needed.

1571:    Developers Note:
1572:    First tries to duplicate the rhs and solution vectors of the `KSP`, if they do not exist tries to get them from the matrix, if
1573:    that does not exist tries to get them from the `DM` (if it is provided).

1575: .seealso: [](chapter_ksp), `MatCreateVecs()`, `VecDestroyVecs()`, `KSPSetWorkVecs()`
1576: @*/
1577: PetscErrorCode KSPCreateVecs(KSP ksp, PetscInt rightn, Vec **right, PetscInt leftn, Vec **left)
1578: {
1579:   Vec       vecr = NULL, vecl = NULL;
1580:   PetscBool matset, pmatset, isshell, preferdm = PETSC_FALSE;
1581:   Mat       mat = NULL;

1583:   if (ksp->dm) {
1584:     PetscObjectTypeCompare((PetscObject)ksp->dm, DMSHELL, &isshell);
1585:     preferdm = isshell ? PETSC_FALSE : PETSC_TRUE;
1586:   }
1587:   if (rightn) {
1589:     if (ksp->vec_sol) vecr = ksp->vec_sol;
1590:     else {
1591:       if (preferdm) {
1592:         DMGetGlobalVector(ksp->dm, &vecr);
1593:       } else if (ksp->pc) {
1594:         PCGetOperatorsSet(ksp->pc, &matset, &pmatset);
1595:         /* check for mat before pmat because for KSPLSQR pmat may be a different size than mat since pmat maybe mat'*mat */
1596:         if (matset) {
1597:           PCGetOperators(ksp->pc, &mat, NULL);
1598:           MatCreateVecs(mat, &vecr, NULL);
1599:         } else if (pmatset) {
1600:           PCGetOperators(ksp->pc, NULL, &mat);
1601:           MatCreateVecs(mat, &vecr, NULL);
1602:         }
1603:       }
1604:       if (!vecr && ksp->dm) DMGetGlobalVector(ksp->dm, &vecr);
1606:     }
1607:     VecDuplicateVecs(vecr, rightn, right);
1608:     if (!ksp->vec_sol) {
1609:       if (preferdm) {
1610:         DMRestoreGlobalVector(ksp->dm, &vecr);
1611:       } else if (mat) {
1612:         VecDestroy(&vecr);
1613:       } else if (ksp->dm) {
1614:         DMRestoreGlobalVector(ksp->dm, &vecr);
1615:       }
1616:     }
1617:   }
1618:   if (leftn) {
1620:     if (ksp->vec_rhs) vecl = ksp->vec_rhs;
1621:     else {
1622:       if (preferdm) {
1623:         DMGetGlobalVector(ksp->dm, &vecl);
1624:       } else if (ksp->pc) {
1625:         PCGetOperatorsSet(ksp->pc, &matset, &pmatset);
1626:         /* check for mat before pmat because for KSPLSQR pmat may be a different size than mat since pmat maybe mat'*mat */
1627:         if (matset) {
1628:           PCGetOperators(ksp->pc, &mat, NULL);
1629:           MatCreateVecs(mat, NULL, &vecl);
1630:         } else if (pmatset) {
1631:           PCGetOperators(ksp->pc, NULL, &mat);
1632:           MatCreateVecs(mat, NULL, &vecl);
1633:         }
1634:       }
1635:       if (!vecl && ksp->dm) DMGetGlobalVector(ksp->dm, &vecl);
1637:     }
1638:     VecDuplicateVecs(vecl, leftn, left);
1639:     if (!ksp->vec_rhs) {
1640:       if (preferdm) {
1641:         DMRestoreGlobalVector(ksp->dm, &vecl);
1642:       } else if (mat) {
1643:         VecDestroy(&vecl);
1644:       } else if (ksp->dm) {
1645:         DMRestoreGlobalVector(ksp->dm, &vecl);
1646:       }
1647:     }
1648:   }
1649:   return 0;
1650: }

1652: /*@C
1653:   KSPSetWorkVecs - Sets a number of work vectors into a `KSP` object

1655:   Collective

1657:   Input Parameters:
1658: + ksp  - iterative context
1659: - nw   - number of work vectors to allocate

1661:   Level: developer

1663:   Developers Note:
1664:   This is PETSC_EXTERN because it may be used by user written plugin KSP implementations

1666: .seealso: [](chapter_ksp), `KSP`, `KSPCreateVecs()`
1667: @*/
1668: PetscErrorCode KSPSetWorkVecs(KSP ksp, PetscInt nw)
1669: {
1670:   VecDestroyVecs(ksp->nwork, &ksp->work);
1671:   ksp->nwork = nw;
1672:   KSPCreateVecs(ksp, nw, &ksp->work, 0, NULL);
1673:   return 0;
1674: }

1676: /*
1677:   KSPDestroyDefault - Destroys a iterative context variable for methods with no separate context.  Preferred calling sequence `KSPDestroy()`.

1679:   Input Parameter:
1680: . ksp - the iterative context

1682:    Developers Note:
1683:    This is PETSC_EXTERN because it may be used by user written plugin KSP implementations

1685: */
1686: PetscErrorCode KSPDestroyDefault(KSP ksp)
1687: {
1689:   PetscFree(ksp->data);
1690:   return 0;
1691: }

1693: /*@
1694:    KSPGetConvergedReason - Gets the reason the `KSP` iteration was stopped.

1696:    Not Collective

1698:    Input Parameter:
1699: .  ksp - the `KSP` context

1701:    Output Parameter:
1702: .  reason - negative value indicates diverged, positive value converged, see KSPConvergedReason

1704:    Possible values for reason: See also manual page for each reason
1705: .vb
1706:    KSP_CONVERGED_RTOL (residual 2-norm decreased by a factor of rtol, from 2-norm of right hand side)
1707:    KSP_CONVERGED_ATOL (residual 2-norm less than abstol)
1708:    KSP_CONVERGED_ITS (used by the preonly preconditioner that always uses ONE iteration, or when the KSPConvergedSkip() convergence test routine is set.
1709:    KSP_CONVERGED_CG_NEG_CURVE (see note below)
1710:    KSP_CONVERGED_CG_CONSTRAINED (see note below)
1711:    KSP_CONVERGED_STEP_LENGTH (see note below)
1712:    KSP_CONVERGED_ITERATING (returned if the solver is not yet finished)
1713:    KSP_DIVERGED_ITS  (required more than its to reach convergence)
1714:    KSP_DIVERGED_DTOL (residual norm increased by a factor of divtol)
1715:    KSP_DIVERGED_NANORINF (residual norm became Not-a-number or Inf likely due to 0/0)
1716:    KSP_DIVERGED_BREAKDOWN (generic breakdown in method)
1717:    KSP_DIVERGED_BREAKDOWN_BICG (Initial residual is orthogonal to preconditioned initial residual. Try a different preconditioner, or a different initial Level.)
1718: .ve

1720:    Options Database Key:
1721: .   -ksp_converged_reason - prints the reason to standard out

1723:    Level: intermediate

1725:    Notes:
1726:     If this routine is called before or doing the `KSPSolve()` the value of `KSP_CONVERGED_ITERATING` is returned

1728:    The values  `KSP_CONVERGED_CG_NEG_CURVE`, `KSP_CONVERGED_CG_CONSTRAINED`, and `KSP_CONVERGED_STEP_LENGTH` are returned only by the special `KSPNASH`, `KSPSTCG`, and `KSPGLTR`
1729:    solvers which are used by the `SNESNEWTONTR` (trust region) solver.

1731: .seealso: [](chapter_ksp), `KSPConvergedReason`, `KSP`, `KSPSetConvergenceTest()`, `KSPConvergedDefault()`, `KSPSetTolerances()`, `KSPConvergedReason`,
1732:           `KSPConvergedReasonView()`
1733: @*/
1734: PetscErrorCode KSPGetConvergedReason(KSP ksp, KSPConvergedReason *reason)
1735: {
1738:   *reason = ksp->reason;
1739:   return 0;
1740: }

1742: /*@C
1743:    KSPGetConvergedReasonString - Return a human readable string for a `KSPConvergedReason`

1745:    Not Collective

1747:    Input Parameter:
1748: .  ksp - the `KSP` context

1750:    Output Parameter:
1751: .  strreason - a human readable string that describes ksp converged reason

1753:    Level: beginner

1755: .seealso: [](chapter_ksp), `KSP`, `KSPGetConvergedReason()`
1756: @*/
1757: PetscErrorCode KSPGetConvergedReasonString(KSP ksp, const char **strreason)
1758: {
1761:   *strreason = KSPConvergedReasons[ksp->reason];
1762:   return 0;
1763: }

1765: #include <petsc/private/dmimpl.h>
1766: /*@
1767:    KSPSetDM - Sets the `DM` that may be used by some preconditioners

1769:    Logically Collective

1771:    Input Parameters:
1772: +  ksp - the preconditioner context
1773: -  dm - the dm, cannot be NULL

1775:    Notes:
1776:    If this is used then the `KSP` will attempt to use the `DM` to create the matrix and use the routine set with
1777:    `DMKSPSetComputeOperators()`. Use `KSPSetDMActive`(ksp,`PETSC_FALSE`) to instead use the matrix you've provided with
1778:    `KSPSetOperators()`.

1780:    A `DM` can only be used for solving one problem at a time because information about the problem is stored on the `DM`,
1781:    even when not using interfaces like `DMKSPSetComputeOperators()`.  Use `DMClone()` to get a distinct `DM` when solving
1782:    different problems using the same function space.

1784:    Level: intermediate

1786: .seealso: [](chapter_ksp), `KSP`, `DM`, `KSPGetDM()`, `KSPSetDMActive()`, `KSPSetComputeOperators()`, `KSPSetComputeRHS()`, `KSPSetComputeInitialGuess()`, `DMKSPSetComputeOperators()`, `DMKSPSetComputeRHS()`, `DMKSPSetComputeInitialGuess()`
1787: @*/
1788: PetscErrorCode KSPSetDM(KSP ksp, DM dm)
1789: {
1790:   PC pc;

1794:   PetscObjectReference((PetscObject)dm);
1795:   if (ksp->dm) { /* Move the DMSNES context over to the new DM unless the new DM already has one */
1796:     if (ksp->dm->dmksp && !dm->dmksp) {
1797:       DMKSP kdm;
1798:       DMCopyDMKSP(ksp->dm, dm);
1799:       DMGetDMKSP(ksp->dm, &kdm);
1800:       if (kdm->originaldm == ksp->dm) kdm->originaldm = dm; /* Grant write privileges to the replacement DM */
1801:     }
1802:     DMDestroy(&ksp->dm);
1803:   }
1804:   ksp->dm     = dm;
1805:   ksp->dmAuto = PETSC_FALSE;
1806:   KSPGetPC(ksp, &pc);
1807:   PCSetDM(pc, dm);
1808:   ksp->dmActive = PETSC_TRUE;
1809:   return 0;
1810: }

1812: /*@
1813:    KSPSetDMActive - Indicates the `DM` should be used to generate the linear system matrix and right hand side vector

1815:    Logically Collective

1817:    Input Parameters:
1818: +  ksp - the preconditioner context
1819: -  flg - use the `DM`

1821:    Note:
1822:    By default `KSPSetDM()` sets the `DM` as active, call `KSPSetDMActive`(ksp,`PETSC_FALSE`); after `KSPSetDM`(ksp,dm) to not have the `KSP` object use the `DM` to generate the matrices.

1824:    Level: intermediate

1826: .seealso: [](chapter_ksp), `KSP`, `DM`, `KSPGetDM()`, `KSPSetDM()`, `SNESSetDM()`, `KSPSetComputeOperators()`, `KSPSetComputeRHS()`, `KSPSetComputeInitialGuess()`
1827: @*/
1828: PetscErrorCode KSPSetDMActive(KSP ksp, PetscBool flg)
1829: {
1832:   ksp->dmActive = flg;
1833:   return 0;
1834: }

1836: /*@
1837:    KSPGetDM - Gets the `DM` that may be used by some preconditioners

1839:    Not Collective

1841:    Input Parameter:
1842: . ksp - the preconditioner context

1844:    Output Parameter:
1845: .  dm - the dm

1847:    Level: intermediate

1849: .seealso: [](chapter_ksp), `KSP`, `DM`, `KSPSetDM()`, `KSPSetDMActive()`
1850: @*/
1851: PetscErrorCode KSPGetDM(KSP ksp, DM *dm)
1852: {
1854:   if (!ksp->dm) {
1855:     DMShellCreate(PetscObjectComm((PetscObject)ksp), &ksp->dm);
1856:     ksp->dmAuto = PETSC_TRUE;
1857:   }
1858:   *dm = ksp->dm;
1859:   return 0;
1860: }

1862: /*@
1863:    KSPSetApplicationContext - Sets the optional user-defined context for the linear solver.

1865:    Logically Collective

1867:    Input Parameters:
1868: +  ksp - the `KSP` context
1869: -  usrP - optional user context

1871:    Fortran Notes:
1872:     To use this from Fortran you must write a Fortran interface definition for this
1873:     function that tells Fortran the Fortran derived data type that you are passing in as the ctx argument.

1875:    Level: intermediate

1877:    Note:
1878:    The user context is a way for users to attach any information to the `KSP` that they may need later when interacting with the `KSP`

1880:    Use `KSPGetApplicationContext()` to get access to the context at a later time.

1882: .seealso: [](chapter_ksp), `KSP`, `KSPGetApplicationContext()`
1883: @*/
1884: PetscErrorCode KSPSetApplicationContext(KSP ksp, void *usrP)
1885: {
1886:   PC pc;

1889:   ksp->user = usrP;
1890:   KSPGetPC(ksp, &pc);
1891:   PCSetApplicationContext(pc, usrP);
1892:   return 0;
1893: }

1895: /*@
1896:    KSPGetApplicationContext - Gets the user-defined context for the linear solver set with `KSPSetApplicationContext()`

1898:    Not Collective

1900:    Input Parameter:
1901: .  ksp - `KSP` context

1903:    Output Parameter:
1904: .  usrP - user context

1906:    Fortran Notes:
1907:     To use this from Fortran you must write a Fortran interface definition for this
1908:     function that tells Fortran the Fortran derived data type that you are passing in as the ctx argument.

1910:    Level: intermediate

1912: .seealso: [](chapter_ksp), `KSP`, `KSPSetApplicationContext()`
1913: @*/
1914: PetscErrorCode KSPGetApplicationContext(KSP ksp, void *usrP)
1915: {
1917:   *(void **)usrP = ksp->user;
1918:   return 0;
1919: }

1921: #include <petsc/private/pcimpl.h>

1923: /*@
1924:    KSPCheckSolve - Checks if the `PCSetUp()` or `KSPSolve()` failed and set the error flag for the outer `PC`. A `KSP_DIVERGED_ITS` is
1925:          not considered a failure in this context

1927:    Collective

1929:    Input Parameters:
1930: +  ksp - the linear solver `KSP` context.
1931: .  pc - the preconditioner context
1932: -  vec - a vector that will be initialized with Inf to indicate lack of convergence

1934:    Note:
1935:    This is called within `KSPSolve()` and `PCApply()` to check if an error has been detected on any particular MPI ranks. By initializing the vector with Inf
1936:    the next call to `KSPCheckNorm()` or `KSPCheckDot()` will provide the same information to all the MPI ranks that an error occurred on at least one of the ranks.

1938:    This may be called by a subset of the processes in the `PC`.

1940:    Level: developer

1942:    Developer Note:
1943:    This is used to manage returning from preconditioners whose inner `KSP` solvers have failed in some way

1945: .seealso: [](chapter_ksp), `KSP`, `KSPCreate()`, `KSPSetType()`, `KSP`, `KSPCheckNorm()`, `KSPCheckDot()`
1946: @*/
1947: PetscErrorCode KSPCheckSolve(KSP ksp, PC pc, Vec vec)
1948: {
1949:   PCFailedReason pcreason;
1950:   PC             subpc;

1952:   KSPGetPC(ksp, &subpc);
1953:   PCGetFailedReason(subpc, &pcreason);
1954:   if (pcreason || (ksp->reason < 0 && ksp->reason != KSP_DIVERGED_ITS)) {
1956:     PetscInfo(ksp, "Detected not converged in KSP inner solve: KSP reason %s PC reason %s\n", KSPConvergedReasons[ksp->reason], PCFailedReasons[pcreason]);
1957:     pc->failedreason = PC_SUBPC_ERROR;
1958:     if (vec) VecSetInf(vec);
1959:   }
1960:   return 0;
1961: }