Actual source code: iterativ.c

petsc-3.14.6 2021-03-30
Report Typos and Errors
  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>

 11: /*@
 12:    KSPGetResidualNorm - Gets the last (approximate preconditioned)
 13:    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: .seealso: KSPBuildResidual()
 26: @*/
 27: PetscErrorCode  KSPGetResidualNorm(KSP ksp,PetscReal *rnorm)
 28: {
 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: {
 60:   *its = ksp->its;
 61:   return(0);
 62: }

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

 67:    Not Collective

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

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

 75:    Level: intermediate

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

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

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

 97:     Logically Collective on ksp

 99:     Input Parameters:
100: +   ksp - the iterative context
101: .   n  - the iteration
102: -   rnorm - the two norm of the residual

104:     Options Database Key:
105: .   -ksp_monitor_singular_value - Activates KSPMonitorSingularValue()

107:     Notes:
108:     The CG solver uses the Lanczos technique for eigenvalue computation,
109:     while GMRES uses the Arnoldi technique; other iterative methods do
110:     not currently compute singular values.

112:     Level: intermediate

114: .seealso: KSPComputeExtremeSingularValues()
115: @*/
116: PetscErrorCode  KSPMonitorSingularValue(KSP ksp,PetscInt n,PetscReal rnorm,PetscViewerAndFormat *dummy)
117: {
118:   PetscReal      emin,emax,c;
120:   PetscViewer    viewer = dummy->viewer;

125:   PetscViewerPushFormat(viewer,dummy->format);
126:   PetscViewerASCIIAddTab(viewer,((PetscObject)ksp)->tablevel);
127:   if (!ksp->calc_sings) {
128:     PetscViewerASCIIPrintf(viewer,"%3D KSP Residual norm %14.12e \n",n,(double)rnorm);
129:   } else {
130:     KSPComputeExtremeSingularValues(ksp,&emax,&emin);
131:     c    = emax/emin;
132:     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)c);
133:   }
134:   PetscViewerASCIISubtractTab(viewer,((PetscObject)ksp)->tablevel);
135:   PetscViewerPopFormat(viewer);
136:   return(0);
137: }

139: /*@C
140:    KSPMonitorSolution - Monitors progress of the KSP solvers by calling
141:    VecView() for the approximate solution at each iteration.

143:    Collective on ksp

145:    Input Parameters:
146: +  ksp - the KSP context
147: .  its - iteration number
148: .  fgnorm - 2-norm of residual (or gradient)
149: -  dummy - a viewer

151:    Level: intermediate

153:    Notes:
154:     For some Krylov methods such as GMRES constructing the solution at
155:   each iteration is expensive, hence using this will slow the code.

157: .seealso: KSPMonitorSet(), KSPMonitorDefault(), VecView()
158: @*/
159: PetscErrorCode  KSPMonitorSolution(KSP ksp,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *dummy)
160: {
162:   Vec            x;
163:   PetscViewer    viewer = dummy->viewer;

167:   KSPBuildSolution(ksp,NULL,&x);
168:   PetscViewerPushFormat(viewer,dummy->format);
169:   VecView(x,viewer);
170:   PetscViewerPopFormat(viewer);
171:   return(0);
172: }

174: /*@C
175:    KSPMonitorDefault - Print the residual norm at each iteration of an
176:    iterative solver.

178:    Collective on ksp

180:    Input Parameters:
181: +  ksp   - iterative context
182: .  n     - iteration number
183: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).
184: -  dummy - an ASCII PetscViewer

186:    Level: intermediate

188: .seealso: KSPMonitorSet(), KSPMonitorTrueResidualNorm(), KSPMonitorLGResidualNormCreate()
189: @*/
190: PetscErrorCode  KSPMonitorDefault(KSP ksp,PetscInt n,PetscReal rnorm,PetscViewerAndFormat *dummy)
191: {
193:   PetscViewer    viewer =  dummy->viewer;

197:   PetscViewerPushFormat(viewer,dummy->format);
198:   PetscViewerASCIIAddTab(viewer,((PetscObject)ksp)->tablevel);
199:   if (n == 0 && ((PetscObject)ksp)->prefix) {
200:     PetscViewerASCIIPrintf(viewer,"  Residual norms for %s solve.\n",((PetscObject)ksp)->prefix);
201:   }
202:   PetscViewerASCIIPrintf(viewer,"%3D KSP Residual norm %14.12e \n",n,(double)rnorm);
203:   PetscViewerASCIISubtractTab(viewer,((PetscObject)ksp)->tablevel);
204:   PetscViewerPopFormat(viewer);
205:   return(0);
206: }

208: /*@C
209:    KSPMonitorTrueResidualNorm - Prints the true residual norm as well as the preconditioned
210:    residual norm at each iteration of an iterative solver.

212:    Collective on ksp

214:    Input Parameters:
215: +  ksp   - iterative context
216: .  n     - iteration number
217: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).
218: -  dummy - an ASCII PetscViewer

220:    Options Database Key:
221: .  -ksp_monitor_true_residual - Activates KSPMonitorTrueResidualNorm()

223:    Notes:
224:    When using right preconditioning, these values are equivalent.

226:    Level: intermediate

228: .seealso: KSPMonitorSet(), KSPMonitorDefault(), KSPMonitorLGResidualNormCreate(),KSPMonitorTrueResidualMaxNorm()
229: @*/
230: PetscErrorCode  KSPMonitorTrueResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,PetscViewerAndFormat *dummy)
231: {
233:   Vec            resid;
234:   PetscReal      truenorm,bnorm;
235:   PetscViewer    viewer = dummy->viewer;
236:   char           normtype[256];

240:   PetscViewerPushFormat(viewer,dummy->format);
241:   PetscViewerASCIIAddTab(viewer,((PetscObject)ksp)->tablevel);
242:   if (n == 0 && ((PetscObject)ksp)->prefix) {
243:     PetscViewerASCIIPrintf(viewer,"  Residual norms for %s solve.\n",((PetscObject)ksp)->prefix);
244:   }
245:   KSPBuildResidual(ksp,NULL,NULL,&resid);
246:   VecNorm(resid,NORM_2,&truenorm);
247:   VecDestroy(&resid);
248:   VecNorm(ksp->vec_rhs,NORM_2,&bnorm);
249:   PetscStrncpy(normtype,KSPNormTypes[ksp->normtype],sizeof(normtype));
250:   PetscStrtolower(normtype);
251:   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));
252:   PetscViewerASCIISubtractTab(viewer,((PetscObject)ksp)->tablevel);
253:   PetscViewerPopFormat(viewer);
254:   return(0);
255: }

257: /*@C
258:    KSPMonitorTrueResidualMaxNorm - Prints the true residual max norm each iteration of an iterative solver.

260:    Collective on ksp

262:    Input Parameters:
263: +  ksp   - iterative context
264: .  n     - iteration number
265: .  rnorm - norm (preconditioned) residual value (may be estimated).
266: -  dummy - an ASCII viewer

268:    Options Database Key:
269: .  -ksp_monitor_max - Activates KSPMonitorTrueResidualMaxNorm()

271:    Notes:
272:    This could be implemented (better) with a flag in ksp.

274:    Level: intermediate

276: .seealso: KSPMonitorSet(), KSPMonitorDefault(), KSPMonitorLGResidualNormCreate(),KSPMonitorTrueResidualNorm()
277: @*/
278: PetscErrorCode  KSPMonitorTrueResidualMaxNorm(KSP ksp,PetscInt n,PetscReal rnorm,PetscViewerAndFormat *dummy)
279: {
281:   Vec            resid;
282:   PetscReal      truenorm,bnorm;
283:   PetscViewer    viewer = dummy->viewer;
284:   char           normtype[256];

288:   PetscViewerPushFormat(viewer,dummy->format);
289:   PetscViewerASCIIAddTab(viewer,((PetscObject)ksp)->tablevel);
290:   if (n == 0 && ((PetscObject)ksp)->prefix) {
291:     PetscViewerASCIIPrintf(viewer,"  Residual norms (max) for %s solve.\n",((PetscObject)ksp)->prefix);
292:   }
293:   KSPBuildResidual(ksp,NULL,NULL,&resid);
294:   VecNorm(resid,NORM_INFINITY,&truenorm);
295:   VecDestroy(&resid);
296:   VecNorm(ksp->vec_rhs,NORM_INFINITY,&bnorm);
297:   PetscStrncpy(normtype,KSPNormTypes[ksp->normtype],sizeof(normtype));
298:   PetscStrtolower(normtype);
299:   PetscViewerASCIIPrintf(viewer,"%3D KSP true resid max norm %14.12e ||r(i)||/||b|| %14.12e\n",n,(double)truenorm,(double)(truenorm/bnorm));
300:   PetscViewerASCIISubtractTab(viewer,((PetscObject)ksp)->tablevel);
301:   PetscViewerPopFormat(viewer);
302:   return(0);
303: }

305: PetscErrorCode  KSPMonitorRange_Private(KSP ksp,PetscInt it,PetscReal *per)
306: {
308:   Vec            resid;
309:   PetscReal      rmax,pwork;
310:   PetscInt       i,n,N;
311:   const PetscScalar *r;

314:   KSPBuildResidual(ksp,NULL,NULL,&resid);
315:   VecNorm(resid,NORM_INFINITY,&rmax);
316:   VecGetLocalSize(resid,&n);
317:   VecGetSize(resid,&N);
318:   VecGetArrayRead(resid,&r);
319:   pwork = 0.0;
320:   for (i=0; i<n; i++) pwork += (PetscAbsScalar(r[i]) > .20*rmax);
321:   VecRestoreArrayRead(resid,&r);
322:   VecDestroy(&resid);
323:   MPIU_Allreduce(&pwork,per,1,MPIU_REAL,MPIU_SUM,PetscObjectComm((PetscObject)ksp));
324:   *per = *per/N;
325:   return(0);
326: }

328: /*@C
329:    KSPMonitorRange - Prints the percentage of residual elements that are more then 10 percent of the maximum value.

331:    Collective on ksp

333:    Input Parameters:
334: +  ksp   - iterative context
335: .  it    - iteration number
336: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).
337: -  dummy - an ASCII viewer

339:    Options Database Key:
340: .  -ksp_monitor_range - Activates KSPMonitorRange()

342:    Level: intermediate

344: .seealso: KSPMonitorSet(), KSPMonitorDefault(), KSPMonitorLGResidualNormCreate()
345: @*/
346: PetscErrorCode  KSPMonitorRange(KSP ksp,PetscInt it,PetscReal rnorm,PetscViewerAndFormat *dummy)
347: {
348:   PetscErrorCode   ierr;
349:   PetscReal        perc,rel;
350:   PetscViewer      viewer = dummy->viewer;
351:   /* should be in a MonitorRangeContext */
352:   static PetscReal prev;

356:   PetscViewerPushFormat(viewer,dummy->format);
357:   PetscViewerASCIIAddTab(viewer,((PetscObject)ksp)->tablevel);
358:   if (!it) prev = rnorm;
359:   if (it == 0 && ((PetscObject)ksp)->prefix) {
360:     PetscViewerASCIIPrintf(viewer,"  Residual norms for %s solve.\n",((PetscObject)ksp)->prefix);
361:   }
362:   KSPMonitorRange_Private(ksp,it,&perc);

364:   rel  = (prev - rnorm)/prev;
365:   prev = rnorm;
366:   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));
367:   PetscViewerASCIISubtractTab(viewer,((PetscObject)ksp)->tablevel);
368:   PetscViewerPopFormat(viewer);
369:   return(0);
370: }

372: /*@C
373:    KSPMonitorDynamicTolerance - Recompute the inner tolerance in every
374:    outer iteration in an adaptive way.

376:    Collective on ksp

378:    Input Parameters:
379: +  ksp   - iterative context
380: .  n     - iteration number (not used)
381: .  fnorm - the current residual norm
382: -  dummy - some context as a C struct. fields:
383:              coef: a scaling coefficient. default 1.0. can be passed through
384:                    -sub_ksp_dynamic_tolerance_param
385:              bnrm: norm of the right-hand side. store it to avoid repeated calculation

387:    Notes:
388:    This may be useful for a flexibly preconditioner Krylov method to
389:    control the accuracy of the inner solves needed to gaurantee the
390:    convergence of the outer iterations.

392:    Level: advanced

394: .seealso: KSPMonitorDynamicToleranceDestroy()
395: @*/
396: PetscErrorCode KSPMonitorDynamicTolerance(KSP ksp,PetscInt its,PetscReal fnorm,void *dummy)
397: {
399:   PC             pc;
400:   PetscReal      outer_rtol, outer_abstol, outer_dtol, inner_rtol;
401:   PetscInt       outer_maxits,nksp,first,i;
402:   KSPDynTolCtx   *scale   = (KSPDynTolCtx*)dummy;
403:   KSP            *subksp = NULL;
404:   KSP            kspinner;
405:   PetscBool      flg;

408:   KSPGetPC(ksp, &pc);

410:   /* compute inner_rtol */
411:   if (scale->bnrm < 0.0) {
412:     Vec b;
413:     KSPGetRhs(ksp, &b);
414:     VecNorm(b, NORM_2, &(scale->bnrm));
415:   }
416:   KSPGetTolerances(ksp, &outer_rtol, &outer_abstol, &outer_dtol, &outer_maxits);
417:   inner_rtol = PetscMin(scale->coef * scale->bnrm * outer_rtol / fnorm, 0.999);
418:   /*PetscPrintf(PETSC_COMM_WORLD, "        Inner rtol = %g\n", (double)inner_rtol);*/

420:   /* if pc is ksp */
421:   PetscObjectTypeCompare((PetscObject)pc,PCKSP,&flg);
422:   if (flg) {
423:     PCKSPGetKSP(pc, &kspinner);
424:     KSPSetTolerances(kspinner, inner_rtol, outer_abstol, outer_dtol, outer_maxits);
425:     return(0);
426:   }

428:   /* if pc is bjacobi */
429:   PetscObjectTypeCompare((PetscObject)pc,PCBJACOBI,&flg);
430:   if (flg) {
431:     PCBJacobiGetSubKSP(pc, &nksp, &first, &subksp);
432:     if (subksp) {
433:       for (i=0; i<nksp; i++) {
434:         KSPSetTolerances(subksp[i], inner_rtol, outer_abstol, outer_dtol, outer_maxits);
435:       }
436:       return(0);
437:     }
438:   }

440:   /* if pc is deflation*/
441:   PetscObjectTypeCompare((PetscObject)pc,PCDEFLATION,&flg);
442:   if (flg) {
443:     PCDeflationGetCoarseKSP(pc,&kspinner);
444:     KSPSetTolerances(kspinner,inner_rtol,outer_abstol,outer_dtol,PETSC_DEFAULT);
445:     return(0);
446:   }

448:   /* todo: dynamic tolerance may apply to other types of pc too */
449:   return(0);
450: }

452: /*
453:   Destroy the dummy context used in KSPMonitorDynamicTolerance()
454: */
455: PetscErrorCode KSPMonitorDynamicToleranceDestroy(void **dummy)
456: {

460:   PetscFree(*dummy);
461:   return(0);
462: }

464: /*
465:   Default (short) KSP Monitor, same as KSPMonitorDefault() except
466:   it prints fewer digits of the residual as the residual gets smaller.
467:   This is because the later digits are meaningless and are often
468:   different on different machines; by using this routine different
469:   machines will usually generate the same output.

471:   Deprecated: Intentionally has no manual page
472: */
473: PetscErrorCode  KSPMonitorDefaultShort(KSP ksp,PetscInt its,PetscReal fnorm,PetscViewerAndFormat *dummy)
474: {
476:   PetscViewer    viewer = dummy->viewer;

480:   PetscViewerPushFormat(viewer,dummy->format);
481:   PetscViewerASCIIAddTab(viewer,((PetscObject)ksp)->tablevel);
482:   if (its == 0 && ((PetscObject)ksp)->prefix) {
483:     PetscViewerASCIIPrintf(viewer,"  Residual norms for %s solve.\n",((PetscObject)ksp)->prefix);
484:   }

486:   if (fnorm > 1.e-9) {
487:     PetscViewerASCIIPrintf(viewer,"%3D KSP Residual norm %g \n",its,(double)fnorm);
488:   } else if (fnorm > 1.e-11) {
489:     PetscViewerASCIIPrintf(viewer,"%3D KSP Residual norm %5.3e \n",its,(double)fnorm);
490:   } else {
491:     PetscViewerASCIIPrintf(viewer,"%3D KSP Residual norm < 1.e-11\n",its);
492:   }
493:   PetscViewerASCIISubtractTab(viewer,((PetscObject)ksp)->tablevel);
494:   PetscViewerPopFormat(viewer);
495:   return(0);
496: }

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

502:    Collective on ksp

504:    Input Parameters:
505: +  ksp   - iterative context
506: .  n     - iteration number
507: .  rnorm - 2-norm residual value (may be estimated)
508: -  dummy - unused convergence context

510:    Returns:
511: .  reason - KSP_CONVERGED_ITERATING, KSP_CONVERGED_ITS

513:    Notes:
514:    This should be used as the convergence test with the option
515:    KSPSetNormType(ksp,KSP_NORM_NONE), since norms of the residual are
516:    not computed. Convergence is then declared after the maximum number
517:    of iterations have been reached. Useful when one is using CG or
518:    BiCGStab as a smoother.

520:    Level: advanced

522: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPSetNormType()
523: @*/
524: PetscErrorCode  KSPConvergedSkip(KSP ksp,PetscInt n,PetscReal rnorm,KSPConvergedReason *reason,void *dummy)
525: {
529:   *reason = KSP_CONVERGED_ITERATING;
530:   if (n >= ksp->max_it) *reason = KSP_CONVERGED_ITS;
531:   return(0);
532: }


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

538:    Note Collective

540:    Output Parameter:
541: .  ctx - convergence context

543:    Level: intermediate

545: .seealso: KSPConvergedDefault(), KSPConvergedDefaultDestroy(), KSPSetConvergenceTest(), KSPSetTolerances(),
546:           KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUIRNorm(), KSPConvergedDefaultSetUMIRNorm(), KSPConvergedDefaultSetConvergedMaxits()
547: @*/
548: PetscErrorCode  KSPConvergedDefaultCreate(void **ctx)
549: {
550:   PetscErrorCode         ierr;
551:   KSPConvergedDefaultCtx *cctx;

554:   PetscNew(&cctx);
555:   *ctx = cctx;
556:   return(0);
557: }

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

564:    Collective on ksp

566:    Input Parameters:
567: .  ksp   - iterative context

569:    Options Database:
570: .   -ksp_converged_use_initial_residual_norm

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

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

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

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


583:    Level: intermediate

585: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUMIRNorm(), KSPConvergedDefaultSetConvergedMaxits()
586: @*/
587: PetscErrorCode  KSPConvergedDefaultSetUIRNorm(KSP ksp)
588: {
589:   KSPConvergedDefaultCtx *ctx = (KSPConvergedDefaultCtx*) ksp->cnvP;

593:   if (ksp->converged != KSPConvergedDefault) return(0);
594:   if (ctx->mininitialrtol) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONGSTATE,"Cannot use KSPConvergedDefaultSetUIRNorm() and KSPConvergedDefaultSetUMIRNorm() together");
595:   ctx->initialrtol = PETSC_TRUE;
596:   return(0);
597: }

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

604:    Collective on ksp

606:    Input Parameters:
607: .  ksp   - iterative context

609:    Options Database:
610: .   -ksp_converged_use_min_initial_residual_norm

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

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

617:    Level: intermediate

619: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUIRNorm(), KSPConvergedDefaultSetConvergedMaxits()
620: @*/
621: PetscErrorCode  KSPConvergedDefaultSetUMIRNorm(KSP ksp)
622: {
623:   KSPConvergedDefaultCtx *ctx = (KSPConvergedDefaultCtx*) ksp->cnvP;

627:   if (ksp->converged != KSPConvergedDefault) return(0);
628:   if (ctx->initialrtol) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONGSTATE,"Cannot use KSPConvergedDefaultSetUIRNorm() and KSPConvergedDefaultSetUMIRNorm() together");
629:   ctx->mininitialrtol = PETSC_TRUE;
630:   return(0);
631: }

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

636:    Collective on ksp

638:    Input Parameters:
639: +  ksp - iterative context
640: -  flg - boolean flag

642:    Options Database:
643: .   -ksp_converged_maxits

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

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

650:    Level: intermediate

652: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUMIRNorm(), KSPConvergedDefaultSetUIRNorm()
653: @*/
654: PetscErrorCode  KSPConvergedDefaultSetConvergedMaxits(KSP ksp, PetscBool flg)
655: {
656:   KSPConvergedDefaultCtx *ctx = (KSPConvergedDefaultCtx*) ksp->cnvP;

661:   if (ksp->converged != KSPConvergedDefault) return(0);
662:   ctx->convmaxits = flg;
663:   return(0);
664: }

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

669:    Collective on ksp

671:    Input Parameters:
672: +  ksp   - iterative context
673: .  n     - iteration number
674: .  rnorm - residual norm (may be estimated, depending on the method may be the preconditioned residual norm)
675: -  ctx - convergence context which must be created by KSPConvergedDefaultCreate()

677:    Output Parameter:
678: +   positive - if the iteration has converged
679: .   negative - if the iteration has diverged
680: -   KSP_CONVERGED_ITERATING - otherwise.

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

688:    where:
689: +     rtol = relative tolerance,
690: .     abstol = absolute tolerance.
691: .     dtol = divergence tolerance,
692: -     rnorm_0 is the two norm of the right hand side (or the preconditioned norm, depending on what was set with
693:           KSPSetNormType(). When initial guess is non-zero you
694:           can call KSPConvergedDefaultSetUIRNorm() to use the norm of (b - A*(initial guess))
695:           as the starting point for relative norm convergence testing, that is as rnorm_0

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

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

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

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

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

707:    Level: intermediate

709: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(),
710:           KSPConvergedDefaultSetUIRNorm(), KSPConvergedDefaultSetUMIRNorm(), KSPConvergedDefaultSetConvergedMaxits(), KSPConvergedDefaultCreate(), KSPConvergedDefaultDestroy()
711: @*/
712: PetscErrorCode  KSPConvergedDefault(KSP ksp,PetscInt n,PetscReal rnorm,KSPConvergedReason *reason,void *ctx)
713: {
714:   PetscErrorCode         ierr;
715:   KSPConvergedDefaultCtx *cctx = (KSPConvergedDefaultCtx*) ctx;
716:   KSPNormType            normtype;

722:   if (!cctx) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_NULL,"Convergence context must have been created with KSPConvergedDefaultCreate()");
723:   *reason = KSP_CONVERGED_ITERATING;

725:   if (cctx->convmaxits && n >= ksp->max_it) {
726:     *reason = KSP_CONVERGED_ITS;
727:     PetscInfo1(ksp,"Linear solver has converged. Maximum number of iterations reached %D\n",n);
728:     return(0);
729:   }
730:   KSPGetNormType(ksp,&normtype);
731:   if (normtype == KSP_NORM_NONE) return(0);

733:   if (!n) {
734:     /* if user gives initial guess need to compute norm of b */
735:     if (!ksp->guess_zero && !cctx->initialrtol) {
736:       PetscReal snorm = 0.0;
737:       if (ksp->normtype == KSP_NORM_UNPRECONDITIONED || ksp->pc_side == PC_RIGHT) {
738:         PetscInfo(ksp,"user has provided nonzero initial guess, computing 2-norm of RHS\n");
739:         VecNorm(ksp->vec_rhs,NORM_2,&snorm);        /*     <- b'*b */
740:       } else {
741:         Vec z;
742:         /* Should avoid allocating the z vector each time but cannot stash it in cctx because if KSPReset() is called the vector size might change */
743:         VecDuplicate(ksp->vec_rhs,&z);
744:         KSP_PCApply(ksp,ksp->vec_rhs,z);
745:         if (ksp->normtype == KSP_NORM_PRECONDITIONED) {
746:           PetscInfo(ksp,"user has provided nonzero initial guess, computing 2-norm of preconditioned RHS\n");
747:           VecNorm(z,NORM_2,&snorm);                 /*    dp <- b'*B'*B*b */
748:         } else if (ksp->normtype == KSP_NORM_NATURAL) {
749:           PetscScalar norm;
750:           PetscInfo(ksp,"user has provided nonzero initial guess, computing natural norm of RHS\n");
751:           VecDot(ksp->vec_rhs,z,&norm);
752:           snorm = PetscSqrtReal(PetscAbsScalar(norm));                            /*    dp <- b'*B*b */
753:         }
754:         VecDestroy(&z);
755:       }
756:       /* handle special case of zero RHS and nonzero guess */
757:       if (!snorm) {
758:         PetscInfo(ksp,"Special case, user has provided nonzero initial guess and zero RHS\n");
759:         snorm = rnorm;
760:       }
761:       if (cctx->mininitialrtol) ksp->rnorm0 = PetscMin(snorm,rnorm);
762:       else ksp->rnorm0 = snorm;
763:     } else {
764:       ksp->rnorm0 = rnorm;
765:     }
766:     ksp->ttol = PetscMax(ksp->rtol*ksp->rnorm0,ksp->abstol);
767:   }

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

771:   if (PetscIsInfOrNanReal(rnorm)) {
772:     PCFailedReason pcreason;
773:     PetscInt       sendbuf,recvbuf;
774:     PCGetFailedReasonRank(ksp->pc,&pcreason);
775:     sendbuf = (PetscInt)pcreason;
776:     MPI_Allreduce(&sendbuf,&recvbuf,1,MPIU_INT,MPIU_MAX,PetscObjectComm((PetscObject)ksp));
777:     if (recvbuf) {
778:       *reason = KSP_DIVERGED_PC_FAILED;
779:       PCSetFailedReason(ksp->pc,(PCFailedReason)recvbuf);
780:       PetscInfo(ksp,"Linear solver pcsetup fails, declaring divergence \n");
781:     } else {
782:       *reason = KSP_DIVERGED_NANORINF;
783:       PetscInfo(ksp,"Linear solver has created a not a number (NaN) as the residual norm, declaring divergence \n");
784:     }
785:   } else if (rnorm <= ksp->ttol) {
786:     if (rnorm < ksp->abstol) {
787:       PetscInfo3(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);
788:       *reason = KSP_CONVERGED_ATOL;
789:     } else {
790:       if (cctx->initialrtol) {
791:         PetscInfo4(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);
792:       } else {
793:         PetscInfo4(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);
794:       }
795:       *reason = KSP_CONVERGED_RTOL;
796:     }
797:   } else if (rnorm >= ksp->divtol*ksp->rnorm0) {
798:     PetscInfo3(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);
799:     *reason = KSP_DIVERGED_DTOL;
800:   }
801:   return(0);
802: }

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

807:    Not Collective

809:    Input Parameters:
810: .  ctx - convergence context

812:    Level: intermediate

814: .seealso: KSPConvergedDefault(), KSPConvergedDefaultCreate(), KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(),
815:           KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUIRNorm(), KSPConvergedDefaultSetUMIRNorm()
816: @*/
817: PetscErrorCode  KSPConvergedDefaultDestroy(void *ctx)
818: {
819:   PetscErrorCode         ierr;
820:   KSPConvergedDefaultCtx *cctx = (KSPConvergedDefaultCtx*) ctx;

823:   VecDestroy(&cctx->work);
824:   PetscFree(ctx);
825:   return(0);
826: }

828: /*
829:    KSPBuildSolutionDefault - Default code to create/move the solution.

831:    Collective on ksp

833:    Input Parameters:
834: +  ksp - iterative context
835: -  v   - pointer to the user's vector

837:    Output Parameter:
838: .  V - pointer to a vector containing the solution

840:    Level: advanced

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

844: .seealso: KSPGetSolution(), KSPBuildResidualDefault()
845: */
846: PetscErrorCode KSPBuildSolutionDefault(KSP ksp,Vec v,Vec *V)
847: {

851:   if (ksp->pc_side == PC_RIGHT) {
852:     if (ksp->pc) {
853:       if (v) {
854:         KSP_PCApply(ksp,ksp->vec_sol,v); *V = v;
855:       } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Not working with right preconditioner");
856:     } else {
857:       if (v) {
858:         VecCopy(ksp->vec_sol,v); *V = v;
859:       } else *V = ksp->vec_sol;
860:     }
861:   } else if (ksp->pc_side == PC_SYMMETRIC) {
862:     if (ksp->pc) {
863:       if (ksp->transpose_solve) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Not working with symmetric preconditioner and transpose solve");
864:       if (v) {
865:         PCApplySymmetricRight(ksp->pc,ksp->vec_sol,v);
866:         *V = v;
867:       } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Not working with symmetric preconditioner");
868:     } else {
869:       if (v) {
870:         VecCopy(ksp->vec_sol,v); *V = v;
871:       } else *V = ksp->vec_sol;
872:     }
873:   } else {
874:     if (v) {
875:       VecCopy(ksp->vec_sol,v); *V = v;
876:     } else *V = ksp->vec_sol;
877:   }
878:   return(0);
879: }

881: /*
882:    KSPBuildResidualDefault - Default code to compute the residual.

884:    Collecive on ksp

886:    Input Parameters:
887: .  ksp - iterative context
888: .  t   - pointer to temporary vector
889: .  v   - pointer to user vector

891:    Output Parameter:
892: .  V - pointer to a vector containing the residual

894:    Level: advanced

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

898: .seealso: KSPBuildSolutionDefault()
899: */
900: PetscErrorCode KSPBuildResidualDefault(KSP ksp,Vec t,Vec v,Vec *V)
901: {
903:   Mat            Amat,Pmat;

906:   if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
907:   PCGetOperators(ksp->pc,&Amat,&Pmat);
908:   KSPBuildSolution(ksp,t,NULL);
909:   KSP_MatMult(ksp,Amat,t,v);
910:   VecAYPX(v,-1.0,ksp->vec_rhs);
911:   *V   = v;
912:   return(0);
913: }

915: /*@C
916:   KSPCreateVecs - Gets a number of work vectors.

918:   Collective on ksp

920:   Input Parameters:
921: + ksp  - iterative context
922: . rightn  - number of right work vectors
923: - leftn   - number of left work vectors to allocate

925:   Output Parameter:
926: +  right - the array of vectors created
927: -  left - the array of left vectors

929:    Note: The right vector has as many elements as the matrix has columns. The left
930:      vector has as many elements as the matrix has rows.

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

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

937:    Level: advanced

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

941: @*/
942: PetscErrorCode KSPCreateVecs(KSP ksp,PetscInt rightn, Vec **right,PetscInt leftn,Vec **left)
943: {
945:   Vec            vecr = NULL,vecl = NULL;
946:   PetscBool      matset,pmatset;
947:   Mat            mat = NULL;

950:   if (rightn) {
951:     if (!right) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"You asked for right vectors but did not pass a pointer to hold them");
952:     if (ksp->vec_sol) vecr = ksp->vec_sol;
953:     else {
954:       if (ksp->pc) {
955:         PCGetOperatorsSet(ksp->pc,&matset,&pmatset);
956:         /* check for mat before pmat because for KSPLSQR pmat may be a different size than mat since pmat maybe mat'*mat */
957:         if (matset) {
958:           PCGetOperators(ksp->pc,&mat,NULL);
959:           MatCreateVecs(mat,&vecr,NULL);
960:         } else if (pmatset) {
961:           PCGetOperators(ksp->pc,NULL,&mat);
962:           MatCreateVecs(mat,&vecr,NULL);
963:         }
964:       }
965:       if (!vecr) {
966:         if (ksp->dm) {
967:           DMGetGlobalVector(ksp->dm,&vecr);
968:         } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONGSTATE,"You requested a vector from a KSP that cannot provide one");
969:       }
970:     }
971:     VecDuplicateVecs(vecr,rightn,right);
972:     if (!ksp->vec_sol) {
973:       if (mat) {
974:         VecDestroy(&vecr);
975:       } else if (ksp->dm) {
976:         DMRestoreGlobalVector(ksp->dm,&vecr);
977:       }
978:     }
979:   }
980:   if (leftn) {
981:     if (!left) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"You asked for left vectors but did not pass a pointer to hold them");
982:     if (ksp->vec_rhs) vecl = ksp->vec_rhs;
983:     else {
984:       if (ksp->pc) {
985:         PCGetOperatorsSet(ksp->pc,&matset,&pmatset);
986:         /* check for mat before pmat because for KSPLSQR pmat may be a different size than mat since pmat maybe mat'*mat */
987:         if (matset) {
988:           PCGetOperators(ksp->pc,&mat,NULL);
989:           MatCreateVecs(mat,NULL,&vecl);
990:         } else if (pmatset) {
991:           PCGetOperators(ksp->pc,NULL,&mat);
992:           MatCreateVecs(mat,NULL,&vecl);
993:         }
994:       }
995:       if (!vecl) {
996:         if (ksp->dm) {
997:           DMGetGlobalVector(ksp->dm,&vecl);
998:         } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONGSTATE,"You requested a vector from a KSP that cannot provide one");
999:       }
1000:     }
1001:     VecDuplicateVecs(vecl,leftn,left);
1002:     if (!ksp->vec_rhs) {
1003:       if (mat) {
1004:         VecDestroy(&vecl);
1005:       } else if (ksp->dm) {
1006:         DMRestoreGlobalVector(ksp->dm,&vecl);
1007:       }
1008:     }
1009:   }
1010:   return(0);
1011: }

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

1016:   Collective on ksp

1018:   Input Parameters:
1019: + ksp  - iterative context
1020: - nw   - number of work vectors to allocate

1022:   Level: developer

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

1026: @*/
1027: PetscErrorCode KSPSetWorkVecs(KSP ksp,PetscInt nw)
1028: {

1032:   VecDestroyVecs(ksp->nwork,&ksp->work);
1033:   ksp->nwork = nw;
1034:   KSPCreateVecs(ksp,nw,&ksp->work,0,NULL);
1035:   PetscLogObjectParents(ksp,nw,ksp->work);
1036:   return(0);
1037: }

1039: /*
1040:   KSPDestroyDefault - Destroys a iterative context variable for methods with
1041:   no separate context.  Preferred calling sequence KSPDestroy().

1043:   Input Parameter:
1044: . ksp - the iterative context

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

1048: */
1049: PetscErrorCode KSPDestroyDefault(KSP ksp)
1050: {

1055:   PetscFree(ksp->data);
1056:   return(0);
1057: }

1059: /*@
1060:    KSPGetConvergedReason - Gets the reason the KSP iteration was stopped.

1062:    Not Collective

1064:    Input Parameter:
1065: .  ksp - the KSP context

1067:    Output Parameter:
1068: .  reason - negative value indicates diverged, positive value converged, see KSPConvergedReason

1070:    Possible values for reason: See also manual page for each reason
1071: $  KSP_CONVERGED_RTOL (residual 2-norm decreased by a factor of rtol, from 2-norm of right hand side)
1072: $  KSP_CONVERGED_ATOL (residual 2-norm less than abstol)
1073: $  KSP_CONVERGED_ITS (used by the preonly preconditioner that always uses ONE iteration, or when the KSPConvergedSkip() convergence test routine is set.
1074: $  KSP_CONVERGED_CG_NEG_CURVE (see note below)
1075: $  KSP_CONVERGED_CG_CONSTRAINED (see note below)
1076: $  KSP_CONVERGED_STEP_LENGTH (see note below)
1077: $  KSP_CONVERGED_ITERATING (returned if the solver is not yet finished)
1078: $  KSP_DIVERGED_ITS  (required more than its to reach convergence)
1079: $  KSP_DIVERGED_DTOL (residual norm increased by a factor of divtol)
1080: $  KSP_DIVERGED_NANORINF (residual norm became Not-a-number or Inf likely due to 0/0)
1081: $  KSP_DIVERGED_BREAKDOWN (generic breakdown in method)
1082: $  KSP_DIVERGED_BREAKDOWN_BICG (Initial residual is orthogonal to preconditioned initial residual. Try a different preconditioner, or a different initial Level.)

1084:    Options Database:
1085: .   -ksp_converged_reason - prints the reason to standard out

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

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

1093:    Level: intermediate

1095: .seealso: KSPSetConvergenceTest(), KSPConvergedDefault(), KSPSetTolerances(), KSPConvergedReason,
1096:           KSPConvergedReasonView()
1097: @*/
1098: PetscErrorCode  KSPGetConvergedReason(KSP ksp,KSPConvergedReason *reason)
1099: {
1103:   *reason = ksp->reason;
1104:   return(0);
1105: }

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

1111:    Logically Collective on ksp

1113:    Input Parameters:
1114: +  ksp - the preconditioner context
1115: -  dm - the dm, cannot be NULL

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

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

1126:    Level: intermediate

1128: .seealso: KSPGetDM(), KSPSetDMActive(), KSPSetComputeOperators(), KSPSetComputeRHS(), KSPSetComputeInitialGuess(), DMKSPSetComputeOperators(), DMKSPSetComputeRHS(), DMKSPSetComputeInitialGuess()
1129: @*/
1130: PetscErrorCode  KSPSetDM(KSP ksp,DM dm)
1131: {
1133:   PC             pc;

1138:   PetscObjectReference((PetscObject)dm);
1139:   if (ksp->dm) {                /* Move the DMSNES context over to the new DM unless the new DM already has one */
1140:     if (ksp->dm->dmksp && !dm->dmksp) {
1141:       DMKSP kdm;
1142:       DMCopyDMKSP(ksp->dm,dm);
1143:       DMGetDMKSP(ksp->dm,&kdm);
1144:       if (kdm->originaldm == ksp->dm) kdm->originaldm = dm; /* Grant write privileges to the replacement DM */
1145:     }
1146:     DMDestroy(&ksp->dm);
1147:   }
1148:   ksp->dm       = dm;
1149:   ksp->dmAuto   = PETSC_FALSE;
1150:   KSPGetPC(ksp,&pc);
1151:   PCSetDM(pc,dm);
1152:   ksp->dmActive = PETSC_TRUE;
1153:   return(0);
1154: }

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

1159:    Logically Collective on ksp

1161:    Input Parameters:
1162: +  ksp - the preconditioner context
1163: -  flg - use the DM

1165:    Notes:
1166:    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.

1168:    Level: intermediate

1170: .seealso: KSPGetDM(), KSPSetDM(), SNESSetDM(), KSPSetComputeOperators(), KSPSetComputeRHS(), KSPSetComputeInitialGuess()
1171: @*/
1172: PetscErrorCode  KSPSetDMActive(KSP ksp,PetscBool flg)
1173: {
1177:   ksp->dmActive = flg;
1178:   return(0);
1179: }

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

1184:    Not Collective

1186:    Input Parameter:
1187: . ksp - the preconditioner context

1189:    Output Parameter:
1190: .  dm - the dm

1192:    Level: intermediate


1195: .seealso: KSPSetDM(), KSPSetDMActive()
1196: @*/
1197: PetscErrorCode  KSPGetDM(KSP ksp,DM *dm)
1198: {

1203:   if (!ksp->dm) {
1204:     DMShellCreate(PetscObjectComm((PetscObject)ksp),&ksp->dm);
1205:     ksp->dmAuto = PETSC_TRUE;
1206:   }
1207:   *dm = ksp->dm;
1208:   return(0);
1209: }

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

1214:    Logically Collective on ksp

1216:    Input Parameters:
1217: +  ksp - the KSP context
1218: -  usrP - optional user context

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

1224:    Level: intermediate

1226: .seealso: KSPGetApplicationContext()
1227: @*/
1228: PetscErrorCode  KSPSetApplicationContext(KSP ksp,void *usrP)
1229: {
1231:   PC             pc;

1235:   ksp->user = usrP;
1236:   KSPGetPC(ksp,&pc);
1237:   PCSetApplicationContext(pc,usrP);
1238:   return(0);
1239: }

1241: /*@
1242:    KSPGetApplicationContext - Gets the user-defined context for the linear solver.

1244:    Not Collective

1246:    Input Parameter:
1247: .  ksp - KSP context

1249:    Output Parameter:
1250: .  usrP - user context

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

1256:    Level: intermediate

1258: .seealso: KSPSetApplicationContext()
1259: @*/
1260: PetscErrorCode  KSPGetApplicationContext(KSP ksp,void *usrP)
1261: {
1264:   *(void**)usrP = ksp->user;
1265:   return(0);
1266: }

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

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

1274:    Collective on ksp

1276:    Input Parameter:
1277: +  ksp - the linear solver (KSP) context.
1278: .  pc - the preconditioner context
1279: -  vec - a vector that will be initialized with Inf to indicate lack of convergence

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

1283:    Level: developer

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

1287: .seealso: KSPCreate(), KSPSetType(), KSP, KSPCheckNorm(), KSPCheckDot()
1288: @*/
1289: PetscErrorCode KSPCheckSolve(KSP ksp,PC pc,Vec vec)
1290: {
1291:   PetscErrorCode     ierr;
1292:   PCFailedReason     pcreason;
1293:   PC                 subpc;

1296:   KSPGetPC(ksp,&subpc);
1297:   PCGetFailedReason(subpc,&pcreason);
1298:   if (pcreason || (ksp->reason < 0 && ksp->reason != KSP_DIVERGED_ITS)) {
1299:     if (pc->erroriffailure) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_NOT_CONVERGED,"Detected not converged in KSP inner solve: KSP reason %s PC reason %s",KSPConvergedReasons[ksp->reason],PCFailedReasons[pcreason]);
1300:     else {
1301:       PetscInfo2(ksp,"Detected not converged in KSP inner solve: KSP reason %s PC reason %s\n",KSPConvergedReasons[ksp->reason],PCFailedReasons[pcreason]);
1302:       pc->failedreason = PC_SUBPC_ERROR;
1303:       if (vec) {
1304:         VecSetInf(vec);
1305:       }
1306:     }
1307:   }
1308:   return(0);
1309: }