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 (approximate preconditioned)
 14:    residual norm that has been computed.

 16:    Not Collective

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

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

 24:    Level: intermediate

 26: .seealso: KSPBuildResidual()
 27: @*/
 28: PetscErrorCode  KSPGetResidualNorm(KSP ksp,PetscReal *rnorm)
 29: {
 32:   *rnorm = ksp->rnorm;
 33:   return 0;
 34: }

 36: /*@
 37:    KSPGetIterationNumber - Gets the current iteration number; if the
 38:          KSPSolve() is complete, returns the number of iterations
 39:          used.

 41:    Not Collective

 43:    Input Parameters:
 44: .  ksp - the iterative context

 46:    Output Parameters:
 47: .  its - number of iterations

 49:    Level: intermediate

 51:    Notes:
 52:       During the ith iteration this returns i-1
 53: .seealso: KSPBuildResidual(), KSPGetResidualNorm(), KSPGetTotalIterations()
 54: @*/
 55: PetscErrorCode  KSPGetIterationNumber(KSP ksp,PetscInt *its)
 56: {
 59:   *its = ksp->its;
 60:   return 0;
 61: }

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

 66:    Not Collective

 68:    Input Parameters:
 69: .  ksp - the iterative context

 71:    Output Parameters:
 72: .  its - total number of iterations

 74:    Level: intermediate

 76:    Notes:
 77:     Use KSPGetIterationNumber() to get the count for the most recent solve only
 78:    If this is called within a linear solve (such as in a KSPMonitor routine) then it does not include iterations within that current solve

 80: .seealso: KSPBuildResidual(), KSPGetResidualNorm(), KSPGetIterationNumber()
 81: @*/
 82: PetscErrorCode  KSPGetTotalIterations(KSP ksp,PetscInt *its)
 83: {
 86:   *its = ksp->totalits;
 87:   return 0;
 88: }

 90: /*@C
 91:   KSPMonitorResidual - Print the preconditioned residual norm at each iteration of an iterative solver.

 93:   Collective on ksp

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

101:   Options Database Key:
102: . -ksp_monitor - Activates KSPMonitorResidual()

104:   Level: intermediate

106: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
107: @*/
108: PetscErrorCode KSPMonitorResidual(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
109: {
110:   PetscViewer       viewer = vf->viewer;
111:   PetscViewerFormat format = vf->format;
112:   PetscInt          tablevel;
113:   const char       *prefix;

116:   PetscObjectGetTabLevel((PetscObject) ksp, &tablevel);
117:   PetscObjectGetOptionsPrefix((PetscObject) ksp, &prefix);
118:   PetscViewerPushFormat(viewer, format);
119:   PetscViewerASCIIAddTab(viewer, tablevel);
120:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
121:   PetscViewerASCIIPrintf(viewer, "%3D KSP Residual norm %14.12e \n", n, (double) rnorm);
122:   PetscViewerASCIISubtractTab(viewer, tablevel);
123:   PetscViewerPopFormat(viewer);
124:   return 0;
125: }

127: /*@C
128:   KSPMonitorResidualDraw - Plots the preconditioned residual at each iteration of an iterative solver.

130:   Collective on ksp

132:   Input Parameters:
133: + ksp   - iterative context
134: . n     - iteration number
135: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
136: - vf    - The viewer context

138:   Options Database Key:
139: . -ksp_monitor draw - Activates KSPMonitorResidualDraw()

141:   Level: intermediate

143: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
144: @*/
145: PetscErrorCode KSPMonitorResidualDraw(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
146: {
147:   PetscViewer       viewer = vf->viewer;
148:   PetscViewerFormat format = vf->format;
149:   Vec               r;

152:   PetscViewerPushFormat(viewer, format);
153:   KSPBuildResidual(ksp, NULL, NULL, &r);
154:   PetscObjectSetName((PetscObject) r, "Residual");
155:   PetscObjectCompose((PetscObject) r, "__Vec_bc_zero__", (PetscObject) ksp);
156:   VecView(r, viewer);
157:   PetscObjectCompose((PetscObject) r, "__Vec_bc_zero__", NULL);
158:   VecDestroy(&r);
159:   PetscViewerPopFormat(viewer);
160:   return 0;
161: }

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

166:   Collective on ksp

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

174:   Options Database Key:
175: . -ksp_monitor draw::draw_lg - Activates KSPMonitorResidualDrawLG()

177:   Level: intermediate

179: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
180: @*/
181: PetscErrorCode KSPMonitorResidualDrawLG(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
182: {
183:   PetscViewer        viewer = vf->viewer;
184:   PetscViewerFormat  format = vf->format;
185:   PetscDrawLG        lg     = vf->lg;
186:   KSPConvergedReason reason;
187:   PetscReal          x, y;

191:   PetscViewerPushFormat(viewer, format);
192:   if (!n) PetscDrawLGReset(lg);
193:   x = (PetscReal) n;
194:   if (rnorm > 0.0) y = PetscLog10Real(rnorm);
195:   else y = -15.0;
196:   PetscDrawLGAddPoint(lg, &x, &y);
197:   KSPGetConvergedReason(ksp, &reason);
198:   if (n <= 20 || !(n % 5) || reason) {
199:     PetscDrawLGDraw(lg);
200:     PetscDrawLGSave(lg);
201:   }
202:   PetscViewerPopFormat(viewer);
203:   return 0;
204: }

206: /*@C
207:   KSPMonitorResidualDrawLGCreate - Creates the plotter for the preconditioned residual.

209:   Collective on ksp

211:   Input Parameters:
212: + viewer - The PetscViewer
213: . format - The viewer format
214: - ctx    - An optional user context

216:   Output Parameter:
217: . vf    - The viewer context

219:   Level: intermediate

221: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
222: @*/
223: PetscErrorCode KSPMonitorResidualDrawLGCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
224: {
225:   PetscViewerAndFormatCreate(viewer, format, vf);
226:   (*vf)->data = ctx;
227:   KSPMonitorLGCreate(PetscObjectComm((PetscObject) viewer), NULL, NULL, "Log Residual Norm", 1, NULL, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &(*vf)->lg);
228:   return 0;
229: }

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

236:   Deprecated: Intentionally has no manual page
237: */
238: PetscErrorCode KSPMonitorResidualShort(KSP ksp, PetscInt its, PetscReal fnorm, PetscViewerAndFormat *vf)
239: {
240:   PetscViewer       viewer = vf->viewer;
241:   PetscViewerFormat format = vf->format;
242:   PetscInt          tablevel;
243:   const char       *prefix;

246:   PetscObjectGetTabLevel((PetscObject) ksp, &tablevel);
247:   PetscObjectGetOptionsPrefix((PetscObject) ksp, &prefix);
248:   PetscViewerPushFormat(viewer, format);
249:   PetscViewerASCIIAddTab(viewer, tablevel);
250:   if (its == 0 && prefix)  PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
251:   if (fnorm > 1.e-9)       PetscViewerASCIIPrintf(viewer, "%3D KSP Residual norm %g \n", its, (double) fnorm);
252:   else if (fnorm > 1.e-11) PetscViewerASCIIPrintf(viewer, "%3D KSP Residual norm %5.3e \n", its, (double) fnorm);
253:   else                     PetscViewerASCIIPrintf(viewer, "%3D KSP Residual norm < 1.e-11\n", its);
254:   PetscViewerASCIISubtractTab(viewer, tablevel);
255:   PetscViewerPopFormat(viewer);
256:   return 0;
257: }

259: PetscErrorCode KSPMonitorRange_Private(KSP ksp, PetscInt it, PetscReal *per)
260: {
261:   Vec                resid;
262:   const PetscScalar *r;
263:   PetscReal          rmax, pwork;
264:   PetscInt           i, n, N;

266:   KSPBuildResidual(ksp, NULL, NULL, &resid);
267:   VecNorm(resid, NORM_INFINITY, &rmax);
268:   VecGetLocalSize(resid, &n);
269:   VecGetSize(resid, &N);
270:   VecGetArrayRead(resid, &r);
271:   pwork = 0.0;
272:   for (i = 0; i < n; ++i) pwork += (PetscAbsScalar(r[i]) > .20*rmax);
273:   VecRestoreArrayRead(resid, &r);
274:   VecDestroy(&resid);
275:   MPIU_Allreduce(&pwork, per, 1, MPIU_REAL, MPIU_SUM, PetscObjectComm((PetscObject) ksp));
276:   *per = *per/N;
277:   return 0;
278: }

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

283:   Collective on ksp

285:   Input Parameters:
286: + ksp   - iterative context
287: . it    - iteration number
288: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
289: - vf    - The viewer context

291:   Options Database Key:
292: . -ksp_monitor_range - Activates KSPMonitorResidualRange()

294:   Level: intermediate

296: .seealso: KSPMonitorSet(), KSPMonitorResidual()
297: @*/
298: PetscErrorCode KSPMonitorResidualRange(KSP ksp, PetscInt it, PetscReal rnorm, PetscViewerAndFormat *vf)
299: {
300:   static PetscReal  prev;
301:   PetscViewer       viewer = vf->viewer;
302:   PetscViewerFormat format = vf->format;
303:   PetscInt          tablevel;
304:   const char       *prefix;
305:   PetscReal         perc, rel;

308:   PetscObjectGetTabLevel((PetscObject) ksp, &tablevel);
309:   PetscObjectGetOptionsPrefix((PetscObject) ksp, &prefix);
310:   PetscViewerPushFormat(viewer, format);
311:   PetscViewerASCIIAddTab(viewer, tablevel);
312:   if (!it) prev = rnorm;
313:   if (it == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
314:   KSPMonitorRange_Private(ksp, it, &perc);
315:   rel  = (prev - rnorm)/prev;
316:   prev = rnorm;
317:   PetscViewerASCIIPrintf(viewer, "%3D 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));
318:   PetscViewerASCIISubtractTab(viewer, tablevel);
319:   PetscViewerPopFormat(viewer);
320:   return 0;
321: }

323: /*@C
324:   KSPMonitorTrueResidual - Prints the true residual norm, as well as the preconditioned residual norm, at each iteration of an iterative solver.

326:   Collective on ksp

328:   Input Parameters:
329: + ksp   - iterative context
330: . n     - iteration number
331: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
332: - vf    - The viewer context

334:   Options Database Key:
335: . -ksp_monitor_true_residual - Activates KSPMonitorTrueResidual()

337:   Notes:
338:   When using right preconditioning, these values are equivalent.

340:   Level: intermediate

342: .seealso: KSPMonitorSet(), KSPMonitorResidual(),KSPMonitorTrueResidualMaxNorm()
343: @*/
344: PetscErrorCode KSPMonitorTrueResidual(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
345: {
346:   PetscViewer       viewer = vf->viewer;
347:   PetscViewerFormat format = vf->format;
348:   Vec               r;
349:   PetscReal         truenorm, bnorm;
350:   char              normtype[256];
351:   PetscInt          tablevel;
352:   const char       *prefix;

355:   PetscObjectGetTabLevel((PetscObject) ksp, &tablevel);
356:   PetscObjectGetOptionsPrefix((PetscObject) ksp, &prefix);
357:   PetscStrncpy(normtype, KSPNormTypes[ksp->normtype], sizeof(normtype));
358:   PetscStrtolower(normtype);
359:   KSPBuildResidual(ksp, NULL, NULL, &r);
360:   VecNorm(r, NORM_2, &truenorm);
361:   VecNorm(ksp->vec_rhs, NORM_2, &bnorm);
362:   VecDestroy(&r);

364:   PetscViewerPushFormat(viewer, format);
365:   PetscViewerASCIIAddTab(viewer, tablevel);
366:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
367:   PetscViewerASCIIPrintf(viewer, "%3D 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));
368:   PetscViewerASCIISubtractTab(viewer, tablevel);
369:   PetscViewerPopFormat(viewer);
370:   return 0;
371: }

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

376:   Collective on ksp

378:   Input Parameters:
379: + ksp   - iterative context
380: . n     - iteration number
381: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
382: - vf    - The viewer context

384:   Options Database Key:
385: . -ksp_monitor_true_residual draw - Activates KSPMonitorResidualDraw()

387:   Level: intermediate

389: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
390: @*/
391: PetscErrorCode KSPMonitorTrueResidualDraw(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
392: {
393:   PetscViewer       viewer = vf->viewer;
394:   PetscViewerFormat format = vf->format;
395:   Vec               r;

398:   PetscViewerPushFormat(viewer, format);
399:   KSPBuildResidual(ksp, NULL, NULL, &r);
400:   PetscObjectSetName((PetscObject) r, "Residual");
401:   PetscObjectCompose((PetscObject) r, "__Vec_bc_zero__", (PetscObject) ksp);
402:   VecView(r, viewer);
403:   PetscObjectCompose((PetscObject) r, "__Vec_bc_zero__", NULL);
404:   VecDestroy(&r);
405:   PetscViewerPopFormat(viewer);
406:   return 0;
407: }

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

412:   Collective on ksp

414:   Input Parameters:
415: + ksp   - iterative context
416: . n     - iteration number
417: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
418: - vf    - The viewer context

420:   Options Database Key:
421: . -ksp_monitor_true_residual draw::draw_lg - Activates KSPMonitorTrueResidualDrawLG()

423:   Level: intermediate

425: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
426: @*/
427: PetscErrorCode KSPMonitorTrueResidualDrawLG(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
428: {
429:   PetscViewer        viewer = vf->viewer;
430:   PetscViewerFormat  format = vf->format;
431:   PetscDrawLG        lg     = vf->lg;
432:   Vec                r;
433:   KSPConvergedReason reason;
434:   PetscReal          truenorm, x[2], y[2];

438:   KSPBuildResidual(ksp, NULL, NULL, &r);
439:   VecNorm(r, NORM_2, &truenorm);
440:   VecDestroy(&r);
441:   PetscViewerPushFormat(viewer, format);
442:   if (!n) PetscDrawLGReset(lg);
443:   x[0] = (PetscReal) n;
444:   if (rnorm > 0.0) y[0] = PetscLog10Real(rnorm);
445:   else y[0] = -15.0;
446:   x[1] = (PetscReal) n;
447:   if (truenorm > 0.0) y[1] = PetscLog10Real(truenorm);
448:   else y[1] = -15.0;
449:   PetscDrawLGAddPoint(lg, x, y);
450:   KSPGetConvergedReason(ksp, &reason);
451:   if (n <= 20 || !(n % 5) || reason) {
452:     PetscDrawLGDraw(lg);
453:     PetscDrawLGSave(lg);
454:   }
455:   PetscViewerPopFormat(viewer);
456:   return 0;
457: }

459: /*@C
460:   KSPMonitorTrueResidualDrawLGCreate - Creates the plotter for the preconditioned residual.

462:   Collective on ksp

464:   Input Parameters:
465: + viewer - The PetscViewer
466: . format - The viewer format
467: - ctx    - An optional user context

469:   Output Parameter:
470: . vf    - The viewer context

472:   Level: intermediate

474: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
475: @*/
476: PetscErrorCode KSPMonitorTrueResidualDrawLGCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
477: {
478:   const char    *names[] = {"preconditioned", "true"};

480:   PetscViewerAndFormatCreate(viewer, format, vf);
481:   (*vf)->data = ctx;
482:   KSPMonitorLGCreate(PetscObjectComm((PetscObject) viewer), NULL, NULL, "Log Residual Norm", 2, names, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &(*vf)->lg);
483:   return 0;
484: }

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

489:   Collective on ksp

491:   Input Parameters:
492: + ksp   - iterative context
493: . n     - iteration number
494: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
495: - vf    - The viewer context

497:   Options Database Key:
498: . -ksp_monitor_true_residual_max - Activates KSPMonitorTrueResidualMax()

500:   Notes:
501:   When using right preconditioning, these values are equivalent.

503:   Level: intermediate

505: .seealso: KSPMonitorSet(), KSPMonitorResidual(),KSPMonitorTrueResidualMaxNorm()
506: @*/
507: PetscErrorCode KSPMonitorTrueResidualMax(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
508: {
509:   PetscViewer       viewer = vf->viewer;
510:   PetscViewerFormat format = vf->format;
511:   Vec               r;
512:   PetscReal         truenorm, bnorm;
513:   char              normtype[256];
514:   PetscInt          tablevel;
515:   const char       *prefix;

518:   PetscObjectGetTabLevel((PetscObject) ksp, &tablevel);
519:   PetscObjectGetOptionsPrefix((PetscObject) ksp, &prefix);
520:   PetscStrncpy(normtype, KSPNormTypes[ksp->normtype], sizeof(normtype));
521:   PetscStrtolower(normtype);
522:   KSPBuildResidual(ksp, NULL, NULL, &r);
523:   VecNorm(r, NORM_INFINITY, &truenorm);
524:   VecNorm(ksp->vec_rhs, NORM_INFINITY, &bnorm);
525:   VecDestroy(&r);

527:   PetscViewerPushFormat(viewer, format);
528:   PetscViewerASCIIAddTab(viewer, tablevel);
529:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
530:   PetscViewerASCIIPrintf(viewer, "%3D KSP %s true resid max norm %14.12e ||r(i)||/||b|| %14.12e\n", n, normtype, (double) truenorm, (double) (truenorm/bnorm));
531:   PetscViewerASCIISubtractTab(viewer, tablevel);
532:   PetscViewerPopFormat(viewer);
533:   return 0;
534: }

536: /*@C
537:   KSPMonitorError - Prints the error norm, as well as the preconditioned residual norm, at each iteration of an iterative solver.

539:   Collective on ksp

541:   Input Parameters:
542: + ksp   - iterative context
543: . n     - iteration number
544: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
545: - vf    - The viewer context

547:   Options Database Key:
548: . -ksp_monitor_error - Activates KSPMonitorError()

550:   Level: intermediate

552: .seealso: KSPMonitorSet(), KSPMonitorResidual(),KSPMonitorTrueResidualMaxNorm()
553: @*/
554: PetscErrorCode KSPMonitorError(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
555: {
556:   PetscViewer       viewer = vf->viewer;
557:   PetscViewerFormat format = vf->format;
558:   DM                dm;
559:   Vec               sol;
560:   PetscReal        *errors;
561:   PetscInt          Nf, f;
562:   PetscInt          tablevel;
563:   const char       *prefix;

566:   PetscObjectGetTabLevel((PetscObject) ksp, &tablevel);
567:   PetscObjectGetOptionsPrefix((PetscObject) ksp, &prefix);
568:   KSPGetDM(ksp, &dm);
569:   DMGetNumFields(dm, &Nf);
570:   DMGetGlobalVector(dm, &sol);
571:   KSPBuildSolution(ksp, sol, NULL);
572:   /* TODO: Make a different monitor that flips sign for SNES, Newton system is A dx = -b, so we need to negate the solution */
573:   VecScale(sol, -1.0);
574:   PetscCalloc1(Nf, &errors);
575:   DMComputeError(dm, sol, errors, NULL);

577:   PetscViewerPushFormat(viewer, format);
578:   PetscViewerASCIIAddTab(viewer, tablevel);
579:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Error norms for %s solve.\n", prefix);
580:   PetscViewerASCIIPrintf(viewer, "%3D KSP Error norm %s", n, Nf > 1 ? "[" : "");
581:   PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
582:   for (f = 0; f < Nf; ++f) {
583:     if (f > 0) PetscViewerASCIIPrintf(viewer, ", ");
584:     PetscViewerASCIIPrintf(viewer, "%14.12e", (double) errors[f]);
585:   }
586:   PetscViewerASCIIPrintf(viewer, "%s resid norm %14.12e\n", Nf > 1 ? "]" : "", (double) rnorm);
587:   PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
588:   PetscViewerASCIISubtractTab(viewer, tablevel);
589:   PetscViewerPopFormat(viewer);
590:   DMRestoreGlobalVector(dm, &sol);
591:   return 0;
592: }

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

597:   Collective on ksp

599:   Input Parameters:
600: + ksp   - iterative context
601: . n     - iteration number
602: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
603: - vf    - The viewer context

605:   Options Database Key:
606: . -ksp_monitor_error draw - Activates KSPMonitorErrorDraw()

608:   Level: intermediate

610: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
611: @*/
612: PetscErrorCode KSPMonitorErrorDraw(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
613: {
614:   PetscViewer       viewer = vf->viewer;
615:   PetscViewerFormat format = vf->format;
616:   DM                dm;
617:   Vec               sol, e;

620:   PetscViewerPushFormat(viewer, format);
621:   KSPGetDM(ksp, &dm);
622:   DMGetGlobalVector(dm, &sol);
623:   KSPBuildSolution(ksp, sol, NULL);
624:   DMComputeError(dm, sol, NULL, &e);
625:   PetscObjectSetName((PetscObject) e, "Error");
626:   PetscObjectCompose((PetscObject) e, "__Vec_bc_zero__", (PetscObject) ksp);
627:   VecView(e, viewer);
628:   PetscObjectCompose((PetscObject) e, "__Vec_bc_zero__", NULL);
629:   VecDestroy(&e);
630:   DMRestoreGlobalVector(dm, &sol);
631:   PetscViewerPopFormat(viewer);
632:   return 0;
633: }

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

638:   Collective on ksp

640:   Input Parameters:
641: + ksp   - iterative context
642: . n     - iteration number
643: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
644: - vf    - The viewer context

646:   Options Database Key:
647: . -ksp_monitor_error draw::draw_lg - Activates KSPMonitorTrueResidualDrawLG()

649:   Level: intermediate

651: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
652: @*/
653: PetscErrorCode KSPMonitorErrorDrawLG(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
654: {
655:   PetscViewer        viewer = vf->viewer;
656:   PetscViewerFormat  format = vf->format;
657:   PetscDrawLG        lg     = vf->lg;
658:   DM                 dm;
659:   Vec                sol;
660:   KSPConvergedReason reason;
661:   PetscReal         *x, *errors;
662:   PetscInt           Nf, f;

666:   KSPGetDM(ksp, &dm);
667:   DMGetNumFields(dm, &Nf);
668:   DMGetGlobalVector(dm, &sol);
669:   KSPBuildSolution(ksp, sol, NULL);
670:   /* TODO: Make a different monitor that flips sign for SNES, Newton system is A dx = -b, so we need to negate the solution */
671:   VecScale(sol, -1.0);
672:   PetscCalloc2(Nf+1, &x, Nf+1, &errors);
673:   DMComputeError(dm, sol, errors, NULL);

675:   PetscViewerPushFormat(viewer, format);
676:   if (!n) PetscDrawLGReset(lg);
677:   for (f = 0; f < Nf; ++f) {
678:     x[f]      = (PetscReal) n;
679:     errors[f] = errors[f] > 0.0 ? PetscLog10Real(errors[f]) : -15.;
680:   }
681:   x[Nf]      = (PetscReal) n;
682:   errors[Nf] = rnorm > 0.0 ? PetscLog10Real(rnorm) : -15.;
683:   PetscDrawLGAddPoint(lg, x, errors);
684:   KSPGetConvergedReason(ksp, &reason);
685:   if (n <= 20 || !(n % 5) || reason) {
686:     PetscDrawLGDraw(lg);
687:     PetscDrawLGSave(lg);
688:   }
689:   PetscViewerPopFormat(viewer);
690:   return 0;
691: }

693: /*@C
694:   KSPMonitorErrorDrawLGCreate - Creates the plotter for the error and preconditioned residual.

696:   Collective on ksp

698:   Input Parameters:
699: + viewer - The PetscViewer
700: . format - The viewer format
701: - ctx    - An optional user context

703:   Output Parameter:
704: . vf    - The viewer context

706:   Level: intermediate

708: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
709: @*/
710: PetscErrorCode KSPMonitorErrorDrawLGCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
711: {
712:   KSP            ksp = (KSP) ctx;
713:   DM             dm;
714:   char         **names;
715:   PetscInt       Nf, f;

717:   KSPGetDM(ksp, &dm);
718:   DMGetNumFields(dm, &Nf);
719:   PetscMalloc1(Nf+1, &names);
720:   for (f = 0; f < Nf; ++f) {
721:     PetscObject disc;
722:     const char *fname;
723:     char        lname[PETSC_MAX_PATH_LEN];

725:     DMGetField(dm, f, NULL, &disc);
726:     PetscObjectGetName(disc, &fname);
727:     PetscStrncpy(lname, fname, PETSC_MAX_PATH_LEN);
728:     PetscStrlcat(lname, " Error", PETSC_MAX_PATH_LEN);
729:     PetscStrallocpy(lname, &names[f]);
730:   }
731:   PetscStrallocpy("residual", &names[Nf]);
732:   PetscViewerAndFormatCreate(viewer, format, vf);
733:   (*vf)->data = ctx;
734:   KSPMonitorLGCreate(PetscObjectComm((PetscObject) viewer), NULL, NULL, "Log Error Norm", Nf+1, (const char **) names, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &(*vf)->lg);
735:   for (f = 0; f <= Nf; ++f) PetscFree(names[f]);
736:   PetscFree(names);
737:   return 0;
738: }

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

743:   Collective on ksp

745:   Input Parameters:
746: + ksp   - iterative context
747: . n     - iteration number
748: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
749: - vf    - The viewer context

751:   Options Database Key:
752: . -ksp_monitor_solution - Activates KSPMonitorSolution()

754:   Level: intermediate

756: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
757: @*/
758: PetscErrorCode KSPMonitorSolution(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
759: {
760:   PetscViewer       viewer = vf->viewer;
761:   PetscViewerFormat format = vf->format;
762:   Vec               x;
763:   PetscReal         snorm;
764:   PetscInt          tablevel;
765:   const char       *prefix;

768:   KSPBuildSolution(ksp, NULL, &x);
769:   VecNorm(x, NORM_2, &snorm);
770:   PetscObjectGetTabLevel((PetscObject) ksp, &tablevel);
771:   PetscObjectGetOptionsPrefix((PetscObject) ksp, &prefix);
772:   PetscViewerPushFormat(viewer, format);
773:   PetscViewerASCIIAddTab(viewer, tablevel);
774:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Solution norms for %s solve.\n", prefix);
775:   PetscViewerASCIIPrintf(viewer, "%3D KSP Solution norm %14.12e \n", n, (double) snorm);
776:   PetscViewerASCIISubtractTab(viewer, tablevel);
777:   PetscViewerPopFormat(viewer);
778:   return 0;
779: }

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

784:   Collective on ksp

786:   Input Parameters:
787: + ksp   - iterative context
788: . n     - iteration number
789: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
790: - vf    - The viewer context

792:   Options Database Key:
793: . -ksp_monitor_solution draw - Activates KSPMonitorSolutionDraw()

795:   Level: intermediate

797: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
798: @*/
799: PetscErrorCode KSPMonitorSolutionDraw(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
800: {
801:   PetscViewer       viewer = vf->viewer;
802:   PetscViewerFormat format = vf->format;
803:   Vec               x;

806:   KSPBuildSolution(ksp, NULL, &x);
807:   PetscViewerPushFormat(viewer, format);
808:   PetscObjectSetName((PetscObject) x, "Solution");
809:   PetscObjectCompose((PetscObject) x, "__Vec_bc_zero__", (PetscObject) ksp);
810:   VecView(x, viewer);
811:   PetscObjectCompose((PetscObject) x, "__Vec_bc_zero__", NULL);
812:   PetscViewerPopFormat(viewer);
813:   return 0;
814: }

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

819:   Collective on ksp

821:   Input Parameters:
822: + ksp   - iterative context
823: . n     - iteration number
824: . rnorm - 2-norm (preconditioned) residual value (may be estimated).
825: - vf    - The viewer context

827:   Options Database Key:
828: . -ksp_monitor_solution draw::draw_lg - Activates KSPMonitorSolutionDrawLG()

830:   Level: intermediate

832: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
833: @*/
834: PetscErrorCode KSPMonitorSolutionDrawLG(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
835: {
836:   PetscViewer        viewer = vf->viewer;
837:   PetscViewerFormat  format = vf->format;
838:   PetscDrawLG        lg     = vf->lg;
839:   Vec                u;
840:   KSPConvergedReason reason;
841:   PetscReal          snorm, x, y;

845:   KSPBuildSolution(ksp, NULL, &u);
846:   VecNorm(u, NORM_2, &snorm);
847:   PetscViewerPushFormat(viewer, format);
848:   if (!n) PetscDrawLGReset(lg);
849:   x = (PetscReal) n;
850:   if (snorm > 0.0) y = PetscLog10Real(snorm);
851:   else y = -15.0;
852:   PetscDrawLGAddPoint(lg, &x, &y);
853:   KSPGetConvergedReason(ksp, &reason);
854:   if (n <= 20 || !(n % 5) || reason) {
855:     PetscDrawLGDraw(lg);
856:     PetscDrawLGSave(lg);
857:   }
858:   PetscViewerPopFormat(viewer);
859:   return 0;
860: }

862: /*@C
863:   KSPMonitorSolutionDrawLGCreate - Creates the plotter for the solution.

865:   Collective on ksp

867:   Input Parameters:
868: + viewer - The PetscViewer
869: . format - The viewer format
870: - ctx    - An optional user context

872:   Output Parameter:
873: . vf    - The viewer context

875:   Level: intermediate

877: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual()
878: @*/
879: PetscErrorCode KSPMonitorSolutionDrawLGCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
880: {
881:   PetscViewerAndFormatCreate(viewer, format, vf);
882:   (*vf)->data = ctx;
883:   KSPMonitorLGCreate(PetscObjectComm((PetscObject) viewer), NULL, NULL, "Log Solution Norm", 1, NULL, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &(*vf)->lg);
884:   return 0;
885: }

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

890:   Logically Collective on ksp

892:   Input Parameters:
893: + ksp   - the iterative context
894: . n     - the iteration
895: . rnorm - the two norm of the residual
896: - vf    - The viewer context

898:   Options Database Key:
899: . -ksp_monitor_singular_value - Activates KSPMonitorSingularValue()

901:   Notes:
902:   The CG solver uses the Lanczos technique for eigenvalue computation,
903:   while GMRES uses the Arnoldi technique; other iterative methods do
904:   not currently compute singular values.

906:   Level: intermediate

908: .seealso: KSPComputeExtremeSingularValues()
909: @*/
910: PetscErrorCode KSPMonitorSingularValue(KSP ksp, PetscInt n, PetscReal rnorm, PetscViewerAndFormat *vf)
911: {
912:   PetscViewer       viewer = vf->viewer;
913:   PetscViewerFormat format = vf->format;
914:   PetscReal         emin, emax;
915:   PetscInt          tablevel;
916:   const char       *prefix;

920:   PetscObjectGetTabLevel((PetscObject) ksp, &tablevel);
921:   PetscObjectGetOptionsPrefix((PetscObject) ksp, &prefix);
922:   PetscViewerPushFormat(viewer, format);
923:   PetscViewerASCIIAddTab(viewer, tablevel);
924:   if (n == 0 && prefix) PetscViewerASCIIPrintf(viewer, "  Residual norms for %s solve.\n", prefix);
925:   if (!ksp->calc_sings) {
926:     PetscViewerASCIIPrintf(viewer, "%3D KSP Residual norm %14.12e \n", n, (double) rnorm);
927:   } else {
928:     KSPComputeExtremeSingularValues(ksp, &emax, &emin);
929:     PetscViewerASCIIPrintf(viewer, "%3D 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));
930:   }
931:   PetscViewerASCIISubtractTab(viewer, tablevel);
932:   PetscViewerPopFormat(viewer);
933:   return 0;
934: }

936: /*@C
937:   KSPMonitorSingularValueCreate - Creates the singular value monitor.

939:   Collective on ksp

941:   Input Parameters:
942: + viewer - The PetscViewer
943: . format - The viewer format
944: - ctx    - An optional user context

946:   Output Parameter:
947: . vf    - The viewer context

949:   Level: intermediate

951: .seealso: KSPMonitorSet(), KSPMonitorSingularValue()
952: @*/
953: PetscErrorCode KSPMonitorSingularValueCreate(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
954: {
955:   KSP            ksp = (KSP) ctx;

957:   PetscViewerAndFormatCreate(viewer, format, vf);
958:   (*vf)->data = ctx;
959:   KSPSetComputeSingularValues(ksp, PETSC_TRUE);
960:   return 0;
961: }

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

966:    Collective on ksp

968:    Input Parameters:
969: +  ksp   - iterative context
970: .  n     - iteration number (not used)
971: .  fnorm - the current residual norm
972: -  dummy - some context as a C struct. fields:
973:              coef: a scaling coefficient. default 1.0. can be passed through
974:                    -sub_ksp_dynamic_tolerance_param
975:              bnrm: norm of the right-hand side. store it to avoid repeated calculation

977:    Notes:
978:    This may be useful for a flexibly preconditioner Krylov method to
979:    control the accuracy of the inner solves needed to guarantee the
980:    convergence of the outer iterations.

982:    Level: advanced

984: .seealso: KSPMonitorDynamicToleranceDestroy()
985: @*/
986: PetscErrorCode KSPMonitorDynamicTolerance(KSP ksp,PetscInt its,PetscReal fnorm,void *dummy)
987: {
988:   PC             pc;
989:   PetscReal      outer_rtol, outer_abstol, outer_dtol, inner_rtol;
990:   PetscInt       outer_maxits,nksp,first,i;
991:   KSPDynTolCtx   *scale   = (KSPDynTolCtx*)dummy;
992:   KSP            *subksp = NULL;
993:   KSP            kspinner;
994:   PetscBool      flg;

996:   KSPGetPC(ksp, &pc);

998:   /* compute inner_rtol */
999:   if (scale->bnrm < 0.0) {
1000:     Vec b;
1001:     KSPGetRhs(ksp, &b);
1002:     VecNorm(b, NORM_2, &(scale->bnrm));
1003:   }
1004:   KSPGetTolerances(ksp, &outer_rtol, &outer_abstol, &outer_dtol, &outer_maxits);
1005:   inner_rtol = PetscMin(scale->coef * scale->bnrm * outer_rtol / fnorm, 0.999);
1006:   /* PetscCall(PetscPrintf(PETSC_COMM_WORLD, "        Inner rtol = %g\n",
1007:      (double)inner_rtol)); */

1009:   /* if pc is ksp */
1010:   PetscObjectTypeCompare((PetscObject)pc,PCKSP,&flg);
1011:   if (flg) {
1012:     PCKSPGetKSP(pc, &kspinner);
1013:     KSPSetTolerances(kspinner, inner_rtol, outer_abstol, outer_dtol, outer_maxits);
1014:     return 0;
1015:   }

1017:   /* if pc is bjacobi */
1018:   PetscObjectTypeCompare((PetscObject)pc,PCBJACOBI,&flg);
1019:   if (flg) {
1020:     PCBJacobiGetSubKSP(pc, &nksp, &first, &subksp);
1021:     if (subksp) {
1022:       for (i=0; i<nksp; i++) {
1023:         KSPSetTolerances(subksp[i], inner_rtol, outer_abstol, outer_dtol, outer_maxits);
1024:       }
1025:       return 0;
1026:     }
1027:   }

1029:   /* if pc is deflation*/
1030:   PetscObjectTypeCompare((PetscObject)pc,PCDEFLATION,&flg);
1031:   if (flg) {
1032:     PCDeflationGetCoarseKSP(pc,&kspinner);
1033:     KSPSetTolerances(kspinner,inner_rtol,outer_abstol,outer_dtol,PETSC_DEFAULT);
1034:     return 0;
1035:   }

1037:   /* todo: dynamic tolerance may apply to other types of pc too */
1038:   return 0;
1039: }

1041: /*
1042:   Destroy the dummy context used in KSPMonitorDynamicTolerance()
1043: */
1044: PetscErrorCode KSPMonitorDynamicToleranceDestroy(void **dummy)
1045: {
1046:   PetscFree(*dummy);
1047:   return 0;
1048: }

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

1054:    Collective on ksp

1056:    Input Parameters:
1057: +  ksp   - iterative context
1058: .  n     - iteration number
1059: .  rnorm - 2-norm residual value (may be estimated)
1060: -  dummy - unused convergence context

1062:    Returns:
1063: .  reason - KSP_CONVERGED_ITERATING, KSP_CONVERGED_ITS

1065:    Notes:
1066:    This should be used as the convergence test with the option
1067:    KSPSetNormType(ksp,KSP_NORM_NONE), since norms of the residual are
1068:    not computed. Convergence is then declared after the maximum number
1069:    of iterations have been reached. Useful when one is using CG or
1070:    BiCGStab as a smoother.

1072:    Level: advanced

1074: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPSetNormType()
1075: @*/
1076: PetscErrorCode  KSPConvergedSkip(KSP ksp,PetscInt n,PetscReal rnorm,KSPConvergedReason *reason,void *dummy)
1077: {
1080:   *reason = KSP_CONVERGED_ITERATING;
1081:   if (n >= ksp->max_it) *reason = KSP_CONVERGED_ITS;
1082:   return 0;
1083: }

1085: /*@C
1086:    KSPConvergedDefaultCreate - Creates and initializes the space used by the KSPConvergedDefault() function context

1088:    Note Collective

1090:    Output Parameter:
1091: .  ctx - convergence context

1093:    Level: intermediate

1095: .seealso: KSPConvergedDefault(), KSPConvergedDefaultDestroy(), KSPSetConvergenceTest(), KSPSetTolerances(),
1096:           KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUIRNorm(), KSPConvergedDefaultSetUMIRNorm(), KSPConvergedDefaultSetConvergedMaxits()
1097: @*/
1098: PetscErrorCode  KSPConvergedDefaultCreate(void **ctx)
1099: {
1100:   KSPConvergedDefaultCtx *cctx;

1102:   PetscNew(&cctx);
1103:   *ctx = cctx;
1104:   return 0;
1105: }

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

1112:    Collective on ksp

1114:    Input Parameters:
1115: .  ksp   - iterative context

1117:    Options Database:
1118: .   -ksp_converged_use_initial_residual_norm <bool> - Use initial residual norm for computing relative convergence

1120:    Notes:
1121:    Use KSPSetTolerances() to alter the defaults for rtol, abstol, dtol.

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

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

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

1130:    Level: intermediate

1132: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUMIRNorm(), KSPConvergedDefaultSetConvergedMaxits()
1133: @*/
1134: PetscErrorCode  KSPConvergedDefaultSetUIRNorm(KSP ksp)
1135: {
1136:   KSPConvergedDefaultCtx *ctx = (KSPConvergedDefaultCtx*) ksp->cnvP;

1139:   if (ksp->converged != KSPConvergedDefault) return 0;
1141:   ctx->initialrtol = PETSC_TRUE;
1142:   return 0;
1143: }

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

1150:    Collective on ksp

1152:    Input Parameters:
1153: .  ksp   - iterative context

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

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

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

1163:    Level: intermediate

1165: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUIRNorm(), KSPConvergedDefaultSetConvergedMaxits()
1166: @*/
1167: PetscErrorCode  KSPConvergedDefaultSetUMIRNorm(KSP ksp)
1168: {
1169:   KSPConvergedDefaultCtx *ctx = (KSPConvergedDefaultCtx*) ksp->cnvP;

1172:   if (ksp->converged != KSPConvergedDefault) return 0;
1174:   ctx->mininitialrtol = PETSC_TRUE;
1175:   return 0;
1176: }

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

1181:    Collective on ksp

1183:    Input Parameters:
1184: +  ksp - iterative context
1185: -  flg - boolean flag

1187:    Options Database:
1188: .   -ksp_converged_maxits <bool> - Declare convergence if the maximum number of iterations is reached

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

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

1195:    Level: intermediate

1197: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUMIRNorm(), KSPConvergedDefaultSetUIRNorm()
1198: @*/
1199: PetscErrorCode  KSPConvergedDefaultSetConvergedMaxits(KSP ksp, PetscBool flg)
1200: {
1201:   KSPConvergedDefaultCtx *ctx = (KSPConvergedDefaultCtx*) ksp->cnvP;

1205:   if (ksp->converged != KSPConvergedDefault) return 0;
1206:   ctx->convmaxits = flg;
1207:   return 0;
1208: }

1210: /*@C
1211:    KSPConvergedDefault - Determines convergence of the linear iterative solvers by default

1213:    Collective on ksp

1215:    Input Parameters:
1216: +  ksp   - iterative context
1217: .  n     - iteration number
1218: .  rnorm - residual norm (may be estimated, depending on the method may be the preconditioned residual norm)
1219: -  ctx - convergence context which must be created by KSPConvergedDefaultCreate()

1221:    Output Parameter:
1222: .  reason - the convergence reason; it is positive if the iteration has converged,
1223:             negative if the iteration has diverged, and KSP_CONVERGED_ITERATING otherwise

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

1231:    where:
1232: +     rtol - relative tolerance,
1233: .     abstol - absolute tolerance.
1234: .     dtol - divergence tolerance,
1235: -     rnorm_0 - the two norm of the right hand side (or the preconditioned norm, depending on what was set with
1236:           KSPSetNormType(). When initial guess is non-zero you
1237:           can call KSPConvergedDefaultSetUIRNorm() to use the norm of (b - A*(initial guess))
1238:           as the starting point for relative norm convergence testing, that is as rnorm_0

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

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

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

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

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

1250:    Level: intermediate

1252: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(),
1253:           KSPConvergedDefaultSetUIRNorm(), KSPConvergedDefaultSetUMIRNorm(), KSPConvergedDefaultSetConvergedMaxits(), KSPConvergedDefaultCreate(), KSPConvergedDefaultDestroy()
1254: @*/
1255: PetscErrorCode  KSPConvergedDefault(KSP ksp,PetscInt n,PetscReal rnorm,KSPConvergedReason *reason,void *ctx)
1256: {
1257:   KSPConvergedDefaultCtx *cctx = (KSPConvergedDefaultCtx*) ctx;
1258:   KSPNormType            normtype;

1264:   *reason = KSP_CONVERGED_ITERATING;

1266:   if (cctx->convmaxits && n >= ksp->max_it) {
1267:     *reason = KSP_CONVERGED_ITS;
1268:     PetscInfo(ksp,"Linear solver has converged. Maximum number of iterations reached %D\n",n);
1269:     return 0;
1270:   }
1271:   KSPGetNormType(ksp,&normtype);
1272:   if (normtype == KSP_NORM_NONE) return 0;

1274:   if (!n) {
1275:     /* if user gives initial guess need to compute norm of b */
1276:     if (!ksp->guess_zero && !cctx->initialrtol) {
1277:       PetscReal snorm = 0.0;
1278:       if (ksp->normtype == KSP_NORM_UNPRECONDITIONED || ksp->pc_side == PC_RIGHT) {
1279:         PetscInfo(ksp,"user has provided nonzero initial guess, computing 2-norm of RHS\n");
1280:         VecNorm(ksp->vec_rhs,NORM_2,&snorm);        /*     <- b'*b */
1281:       } else {
1282:         Vec z;
1283:         /* Should avoid allocating the z vector each time but cannot stash it in cctx because if KSPReset() is called the vector size might change */
1284:         VecDuplicate(ksp->vec_rhs,&z);
1285:         KSP_PCApply(ksp,ksp->vec_rhs,z);
1286:         if (ksp->normtype == KSP_NORM_PRECONDITIONED) {
1287:           PetscInfo(ksp,"user has provided nonzero initial guess, computing 2-norm of preconditioned RHS\n");
1288:           VecNorm(z,NORM_2,&snorm);                 /*    dp <- b'*B'*B*b */
1289:         } else if (ksp->normtype == KSP_NORM_NATURAL) {
1290:           PetscScalar norm;
1291:           PetscInfo(ksp,"user has provided nonzero initial guess, computing natural norm of RHS\n");
1292:           VecDot(ksp->vec_rhs,z,&norm);
1293:           snorm = PetscSqrtReal(PetscAbsScalar(norm));                            /*    dp <- b'*B*b */
1294:         }
1295:         VecDestroy(&z);
1296:       }
1297:       /* handle special case of zero RHS and nonzero guess */
1298:       if (!snorm) {
1299:         PetscInfo(ksp,"Special case, user has provided nonzero initial guess and zero RHS\n");
1300:         snorm = rnorm;
1301:       }
1302:       if (cctx->mininitialrtol) ksp->rnorm0 = PetscMin(snorm,rnorm);
1303:       else ksp->rnorm0 = snorm;
1304:     } else {
1305:       ksp->rnorm0 = rnorm;
1306:     }
1307:     ksp->ttol = PetscMax(ksp->rtol*ksp->rnorm0,ksp->abstol);
1308:   }

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

1312:   if (PetscIsInfOrNanReal(rnorm)) {
1313:     PCFailedReason pcreason;
1314:     PetscInt       sendbuf,recvbuf;
1315:     PCGetFailedReasonRank(ksp->pc,&pcreason);
1316:     sendbuf = (PetscInt)pcreason;
1317:     MPI_Allreduce(&sendbuf,&recvbuf,1,MPIU_INT,MPIU_MAX,PetscObjectComm((PetscObject)ksp));
1318:     if (recvbuf) {
1319:       *reason = KSP_DIVERGED_PC_FAILED;
1320:       PCSetFailedReason(ksp->pc,(PCFailedReason)recvbuf);
1321:       PetscInfo(ksp,"Linear solver pcsetup fails, declaring divergence \n");
1322:     } else {
1323:       *reason = KSP_DIVERGED_NANORINF;
1324:       PetscInfo(ksp,"Linear solver has created a not a number (NaN) as the residual norm, declaring divergence \n");
1325:     }
1326:   } else if (rnorm <= ksp->ttol) {
1327:     if (rnorm < ksp->abstol) {
1328:       PetscInfo(ksp,"Linear solver has converged. Residual norm %14.12e is less than absolute tolerance %14.12e at iteration %D\n",(double)rnorm,(double)ksp->abstol,n);
1329:       *reason = KSP_CONVERGED_ATOL;
1330:     } else {
1331:       if (cctx->initialrtol) {
1332:         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 %D\n",(double)rnorm,(double)ksp->rtol,(double)ksp->rnorm0,n);
1333:       } else {
1334:         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 %D\n",(double)rnorm,(double)ksp->rtol,(double)ksp->rnorm0,n);
1335:       }
1336:       *reason = KSP_CONVERGED_RTOL;
1337:     }
1338:   } else if (rnorm >= ksp->divtol*ksp->rnorm0) {
1339:     PetscInfo(ksp,"Linear solver is diverging. Initial right hand size norm %14.12e, current residual norm %14.12e at iteration %D\n",(double)ksp->rnorm0,(double)rnorm,n);
1340:     *reason = KSP_DIVERGED_DTOL;
1341:   }
1342:   return 0;
1343: }

1345: /*@C
1346:    KSPConvergedDefaultDestroy - Frees the space used by the KSPConvergedDefault() function context

1348:    Not Collective

1350:    Input Parameters:
1351: .  ctx - convergence context

1353:    Level: intermediate

1355: .seealso: KSPConvergedDefault(), KSPConvergedDefaultCreate(), KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(),
1356:           KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUIRNorm(), KSPConvergedDefaultSetUMIRNorm()
1357: @*/
1358: PetscErrorCode  KSPConvergedDefaultDestroy(void *ctx)
1359: {
1360:   KSPConvergedDefaultCtx *cctx = (KSPConvergedDefaultCtx*) ctx;

1362:   VecDestroy(&cctx->work);
1363:   PetscFree(ctx);
1364:   return 0;
1365: }

1367: /*
1368:    KSPBuildSolutionDefault - Default code to create/move the solution.

1370:    Collective on ksp

1372:    Input Parameters:
1373: +  ksp - iterative context
1374: -  v   - pointer to the user's vector

1376:    Output Parameter:
1377: .  V - pointer to a vector containing the solution

1379:    Level: advanced

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

1383: .seealso: KSPGetSolution(), KSPBuildResidualDefault()
1384: */
1385: PetscErrorCode KSPBuildSolutionDefault(KSP ksp,Vec v,Vec *V)
1386: {
1387:   if (ksp->pc_side == PC_RIGHT) {
1388:     if (ksp->pc) {
1389:       if (v) {
1390:         KSP_PCApply(ksp,ksp->vec_sol,v); *V = v;
1391:       } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Not working with right preconditioner");
1392:     } else {
1393:       if (v) {
1394:         VecCopy(ksp->vec_sol,v); *V = v;
1395:       } else *V = ksp->vec_sol;
1396:     }
1397:   } else if (ksp->pc_side == PC_SYMMETRIC) {
1398:     if (ksp->pc) {
1400:       if (v) {
1401:         PCApplySymmetricRight(ksp->pc,ksp->vec_sol,v);
1402:         *V = v;
1403:       } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Not working with symmetric preconditioner");
1404:     } else {
1405:       if (v) {
1406:         VecCopy(ksp->vec_sol,v); *V = v;
1407:       } else *V = ksp->vec_sol;
1408:     }
1409:   } else {
1410:     if (v) {
1411:       VecCopy(ksp->vec_sol,v); *V = v;
1412:     } else *V = ksp->vec_sol;
1413:   }
1414:   return 0;
1415: }

1417: /*
1418:    KSPBuildResidualDefault - Default code to compute the residual.

1420:    Collecive on ksp

1422:    Input Parameters:
1423: .  ksp - iterative context
1424: .  t   - pointer to temporary vector
1425: .  v   - pointer to user vector

1427:    Output Parameter:
1428: .  V - pointer to a vector containing the residual

1430:    Level: advanced

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

1434: .seealso: KSPBuildSolutionDefault()
1435: */
1436: PetscErrorCode KSPBuildResidualDefault(KSP ksp,Vec t,Vec v,Vec *V)
1437: {
1438:   Mat            Amat,Pmat;

1440:   if (!ksp->pc) KSPGetPC(ksp,&ksp->pc);
1441:   PCGetOperators(ksp->pc,&Amat,&Pmat);
1442:   KSPBuildSolution(ksp,t,NULL);
1443:   KSP_MatMult(ksp,Amat,t,v);
1444:   VecAYPX(v,-1.0,ksp->vec_rhs);
1445:   *V   = v;
1446:   return 0;
1447: }

1449: /*@C
1450:   KSPCreateVecs - Gets a number of work vectors.

1452:   Collective on ksp

1454:   Input Parameters:
1455: + ksp  - iterative context
1456: . rightn  - number of right work vectors
1457: - leftn   - number of left work vectors to allocate

1459:   Output Parameters:
1460: +  right - the array of vectors created
1461: -  left - the array of left vectors

1463:    Note: The right vector has as many elements as the matrix has columns. The left
1464:      vector has as many elements as the matrix has rows.

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

1468:    Developers Note: 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
1469:                     that does not exist tries to get them from the DM (if it is provided).

1471:    Level: advanced

1473: .seealso:   MatCreateVecs(), VecDestroyVecs()

1475: @*/
1476: PetscErrorCode KSPCreateVecs(KSP ksp,PetscInt rightn, Vec **right,PetscInt leftn,Vec **left)
1477: {
1478:   Vec            vecr = NULL,vecl = NULL;
1479:   PetscBool      matset,pmatset,isshell,preferdm = PETSC_FALSE;
1480:   Mat            mat = NULL;

1482:   if (ksp->dm) {
1483:     PetscObjectTypeCompare((PetscObject) ksp->dm, DMSHELL, &isshell);
1484:     preferdm = isshell ? PETSC_FALSE : PETSC_TRUE;
1485:   }
1486:   if (rightn) {
1488:     if (ksp->vec_sol) vecr = ksp->vec_sol;
1489:     else {
1490:       if (preferdm) {
1491:         DMGetGlobalVector(ksp->dm,&vecr);
1492:       } else if (ksp->pc) {
1493:         PCGetOperatorsSet(ksp->pc,&matset,&pmatset);
1494:         /* check for mat before pmat because for KSPLSQR pmat may be a different size than mat since pmat maybe mat'*mat */
1495:         if (matset) {
1496:           PCGetOperators(ksp->pc,&mat,NULL);
1497:           MatCreateVecs(mat,&vecr,NULL);
1498:         } else if (pmatset) {
1499:           PCGetOperators(ksp->pc,NULL,&mat);
1500:           MatCreateVecs(mat,&vecr,NULL);
1501:         }
1502:       }
1503:       if (!vecr && ksp->dm) {
1504:         DMGetGlobalVector(ksp->dm,&vecr);
1505:       }
1507:     }
1508:     VecDuplicateVecs(vecr,rightn,right);
1509:     if (!ksp->vec_sol) {
1510:       if (preferdm) {
1511:         DMRestoreGlobalVector(ksp->dm,&vecr);
1512:       } else if (mat) {
1513:         VecDestroy(&vecr);
1514:       } else if (ksp->dm) {
1515:         DMRestoreGlobalVector(ksp->dm,&vecr);
1516:       }
1517:     }
1518:   }
1519:   if (leftn) {
1521:     if (ksp->vec_rhs) vecl = ksp->vec_rhs;
1522:     else {
1523:       if (preferdm) {
1524:         DMGetGlobalVector(ksp->dm,&vecl);
1525:       } else if (ksp->pc) {
1526:         PCGetOperatorsSet(ksp->pc,&matset,&pmatset);
1527:         /* check for mat before pmat because for KSPLSQR pmat may be a different size than mat since pmat maybe mat'*mat */
1528:         if (matset) {
1529:           PCGetOperators(ksp->pc,&mat,NULL);
1530:           MatCreateVecs(mat,NULL,&vecl);
1531:         } else if (pmatset) {
1532:           PCGetOperators(ksp->pc,NULL,&mat);
1533:           MatCreateVecs(mat,NULL,&vecl);
1534:         }
1535:       }
1536:       if (!vecl && ksp->dm) DMGetGlobalVector(ksp->dm,&vecl);
1538:     }
1539:     VecDuplicateVecs(vecl,leftn,left);
1540:     if (!ksp->vec_rhs) {
1541:       if (preferdm) {
1542:         DMRestoreGlobalVector(ksp->dm,&vecl);
1543:       } else if (mat) {
1544:         VecDestroy(&vecl);
1545:       } else if (ksp->dm) {
1546:         DMRestoreGlobalVector(ksp->dm,&vecl);
1547:       }
1548:     }
1549:   }
1550:   return 0;
1551: }

1553: /*@C
1554:   KSPSetWorkVecs - Sets a number of work vectors into a KSP object

1556:   Collective on ksp

1558:   Input Parameters:
1559: + ksp  - iterative context
1560: - nw   - number of work vectors to allocate

1562:   Level: developer

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

1566: @*/
1567: PetscErrorCode KSPSetWorkVecs(KSP ksp,PetscInt nw)
1568: {
1569:   VecDestroyVecs(ksp->nwork,&ksp->work);
1570:   ksp->nwork = nw;
1571:   KSPCreateVecs(ksp,nw,&ksp->work,0,NULL);
1572:   PetscLogObjectParents(ksp,nw,ksp->work);
1573:   return 0;
1574: }

1576: /*
1577:   KSPDestroyDefault - Destroys a iterative context variable for methods with
1578:   no separate context.  Preferred calling sequence KSPDestroy().

1580:   Input Parameter:
1581: . ksp - the iterative context

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

1585: */
1586: PetscErrorCode KSPDestroyDefault(KSP ksp)
1587: {
1589:   PetscFree(ksp->data);
1590:   return 0;
1591: }

1593: /*@
1594:    KSPGetConvergedReason - Gets the reason the KSP iteration was stopped.

1596:    Not Collective

1598:    Input Parameter:
1599: .  ksp - the KSP context

1601:    Output Parameter:
1602: .  reason - negative value indicates diverged, positive value converged, see KSPConvergedReason

1604:    Possible values for reason: See also manual page for each reason
1605: $  KSP_CONVERGED_RTOL (residual 2-norm decreased by a factor of rtol, from 2-norm of right hand side)
1606: $  KSP_CONVERGED_ATOL (residual 2-norm less than abstol)
1607: $  KSP_CONVERGED_ITS (used by the preonly preconditioner that always uses ONE iteration, or when the KSPConvergedSkip() convergence test routine is set.
1608: $  KSP_CONVERGED_CG_NEG_CURVE (see note below)
1609: $  KSP_CONVERGED_CG_CONSTRAINED (see note below)
1610: $  KSP_CONVERGED_STEP_LENGTH (see note below)
1611: $  KSP_CONVERGED_ITERATING (returned if the solver is not yet finished)
1612: $  KSP_DIVERGED_ITS  (required more than its to reach convergence)
1613: $  KSP_DIVERGED_DTOL (residual norm increased by a factor of divtol)
1614: $  KSP_DIVERGED_NANORINF (residual norm became Not-a-number or Inf likely due to 0/0)
1615: $  KSP_DIVERGED_BREAKDOWN (generic breakdown in method)
1616: $  KSP_DIVERGED_BREAKDOWN_BICG (Initial residual is orthogonal to preconditioned initial residual. Try a different preconditioner, or a different initial Level.)

1618:    Options Database:
1619: .   -ksp_converged_reason - prints the reason to standard out

1621:    Notes:
1622:     If this routine is called before or doing the KSPSolve() the value of KSP_CONVERGED_ITERATING is returned

1624:    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
1625:    solvers which are used by the SNESNEWTONTR (trust region) solver.

1627:    Level: intermediate

1629: .seealso: KSPSetConvergenceTest(), KSPConvergedDefault(), KSPSetTolerances(), KSPConvergedReason,
1630:           KSPConvergedReasonView()
1631: @*/
1632: PetscErrorCode  KSPGetConvergedReason(KSP ksp,KSPConvergedReason *reason)
1633: {
1636:   *reason = ksp->reason;
1637:   return 0;
1638: }

1640: /*@C
1641:    KSPGetConvergedReasonString - Return a human readable string for ksp converged reason

1643:    Not Collective

1645:    Input Parameter:
1646: .  ksp - the KSP context

1648:    Output Parameter:
1649: .  strreason - a human readable string that describes ksp converged reason

1651:    Level: beginner

1653: .seealso: KSPGetConvergedReason()
1654: @*/
1655: PetscErrorCode KSPGetConvergedReasonString(KSP ksp,const char** strreason)
1656: {
1659:   *strreason = KSPConvergedReasons[ksp->reason];
1660:   return 0;
1661: }

1663: #include <petsc/private/dmimpl.h>
1664: /*@
1665:    KSPSetDM - Sets the DM that may be used by some preconditioners

1667:    Logically Collective on ksp

1669:    Input Parameters:
1670: +  ksp - the preconditioner context
1671: -  dm - the dm, cannot be NULL

1673:    Notes:
1674:    If this is used then the KSP will attempt to use the DM to create the matrix and use the routine set with
1675:    DMKSPSetComputeOperators(). Use KSPSetDMActive(ksp,PETSC_FALSE) to instead use the matrix you've provided with
1676:    KSPSetOperators().

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

1682:    Level: intermediate

1684: .seealso: KSPGetDM(), KSPSetDMActive(), KSPSetComputeOperators(), KSPSetComputeRHS(), KSPSetComputeInitialGuess(), DMKSPSetComputeOperators(), DMKSPSetComputeRHS(), DMKSPSetComputeInitialGuess()
1685: @*/
1686: PetscErrorCode  KSPSetDM(KSP ksp,DM dm)
1687: {
1688:   PC             pc;

1692:   PetscObjectReference((PetscObject)dm);
1693:   if (ksp->dm) {                /* Move the DMSNES context over to the new DM unless the new DM already has one */
1694:     if (ksp->dm->dmksp && !dm->dmksp) {
1695:       DMKSP kdm;
1696:       DMCopyDMKSP(ksp->dm,dm);
1697:       DMGetDMKSP(ksp->dm,&kdm);
1698:       if (kdm->originaldm == ksp->dm) kdm->originaldm = dm; /* Grant write privileges to the replacement DM */
1699:     }
1700:     DMDestroy(&ksp->dm);
1701:   }
1702:   ksp->dm       = dm;
1703:   ksp->dmAuto   = PETSC_FALSE;
1704:   KSPGetPC(ksp,&pc);
1705:   PCSetDM(pc,dm);
1706:   ksp->dmActive = PETSC_TRUE;
1707:   return 0;
1708: }

1710: /*@
1711:    KSPSetDMActive - Indicates the DM should be used to generate the linear system matrix and right hand side

1713:    Logically Collective on ksp

1715:    Input Parameters:
1716: +  ksp - the preconditioner context
1717: -  flg - use the DM

1719:    Notes:
1720:    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.

1722:    Level: intermediate

1724: .seealso: KSPGetDM(), KSPSetDM(), SNESSetDM(), KSPSetComputeOperators(), KSPSetComputeRHS(), KSPSetComputeInitialGuess()
1725: @*/
1726: PetscErrorCode  KSPSetDMActive(KSP ksp,PetscBool flg)
1727: {
1730:   ksp->dmActive = flg;
1731:   return 0;
1732: }

1734: /*@
1735:    KSPGetDM - Gets the DM that may be used by some preconditioners

1737:    Not Collective

1739:    Input Parameter:
1740: . ksp - the preconditioner context

1742:    Output Parameter:
1743: .  dm - the dm

1745:    Level: intermediate

1747: .seealso: KSPSetDM(), KSPSetDMActive()
1748: @*/
1749: PetscErrorCode  KSPGetDM(KSP ksp,DM *dm)
1750: {
1752:   if (!ksp->dm) {
1753:     DMShellCreate(PetscObjectComm((PetscObject)ksp),&ksp->dm);
1754:     ksp->dmAuto = PETSC_TRUE;
1755:   }
1756:   *dm = ksp->dm;
1757:   return 0;
1758: }

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

1763:    Logically Collective on ksp

1765:    Input Parameters:
1766: +  ksp - the KSP context
1767: -  usrP - optional user context

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

1773:    Level: intermediate

1775: .seealso: KSPGetApplicationContext()
1776: @*/
1777: PetscErrorCode  KSPSetApplicationContext(KSP ksp,void *usrP)
1778: {
1779:   PC             pc;

1782:   ksp->user = usrP;
1783:   KSPGetPC(ksp,&pc);
1784:   PCSetApplicationContext(pc,usrP);
1785:   return 0;
1786: }

1788: /*@
1789:    KSPGetApplicationContext - Gets the user-defined context for the linear solver.

1791:    Not Collective

1793:    Input Parameter:
1794: .  ksp - KSP context

1796:    Output Parameter:
1797: .  usrP - user context

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

1803:    Level: intermediate

1805: .seealso: KSPSetApplicationContext()
1806: @*/
1807: PetscErrorCode  KSPGetApplicationContext(KSP ksp,void *usrP)
1808: {
1810:   *(void**)usrP = ksp->user;
1811:   return 0;
1812: }

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

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

1820:    Collective on ksp

1822:    Input Parameters:
1823: +  ksp - the linear solver (KSP) context.
1824: .  pc - the preconditioner context
1825: -  vec - a vector that will be initialized with Inf to indicate lack of convergence

1827:    Notes: this may be called by a subset of the processes in the PC

1829:    Level: developer

1831:    Developer Note: this is used to manage returning from preconditioners whose inner KSP solvers have failed in some way

1833: .seealso: KSPCreate(), KSPSetType(), KSP, KSPCheckNorm(), KSPCheckDot()
1834: @*/
1835: PetscErrorCode KSPCheckSolve(KSP ksp,PC pc,Vec vec)
1836: {
1837:   PCFailedReason pcreason;
1838:   PC             subpc;

1840:   KSPGetPC(ksp,&subpc);
1841:   PCGetFailedReason(subpc,&pcreason);
1842:   if (pcreason || (ksp->reason < 0 && ksp->reason != KSP_DIVERGED_ITS)) {
1844:     else {
1845:       PetscInfo(ksp,"Detected not converged in KSP inner solve: KSP reason %s PC reason %s\n",KSPConvergedReasons[ksp->reason],PCFailedReasons[pcreason]);
1846:       pc->failedreason = PC_SUBPC_ERROR;
1847:       if (vec) VecSetInf(vec);
1848:     }
1849:   }
1850:   return 0;
1851: }