Actual source code: iterativ.c

petsc-3.9.3 2018-07-02
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 truely 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: .keywords: KSP, get, residual norm

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

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

 43:    Not Collective

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

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

 51:    Level: intermediate

 53:    Notes:
 54:       During the ith iteration this returns i-1
 55: .keywords: KSP, get, residual norm

 57: .seealso: KSPBuildResidual(), KSPGetResidualNorm(), KSPGetTotalIterations()
 58: @*/
 59: PetscErrorCode  KSPGetIterationNumber(KSP ksp,PetscInt *its)
 60: {
 64:   *its = ksp->its;
 65:   return(0);
 66: }

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

 71:    Not Collective

 73:    Input Parameters:
 74: .  ksp - the iterative context

 76:    Output Parameters:
 77: .  its - total number of iterations

 79:    Level: intermediate

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

 84: .keywords: KSP, get, residual norm

 86: .seealso: KSPBuildResidual(), KSPGetResidualNorm(), KSPGetIterationNumber()
 87: @*/
 88: PetscErrorCode  KSPGetTotalIterations(KSP ksp,PetscInt *its)
 89: {
 93:   *its = ksp->totalits;
 94:   return(0);
 95: }

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

102:     Logically Collective on KSP

104:     Input Parameters:
105: +   ksp - the iterative context
106: .   n  - the iteration
107: -   rnorm - the two norm of the residual

109:     Options Database Key:
110: .   -ksp_monitor_singular_value - Activates KSPMonitorSingularValue()

112:     Notes:
113:     The CG solver uses the Lanczos technique for eigenvalue computation,
114:     while GMRES uses the Arnoldi technique; other iterative methods do
115:     not currently compute singular values.

117:     Level: intermediate

119: .keywords: KSP, CG, default, monitor, extreme, singular values, Lanczos, Arnoldi

121: .seealso: KSPComputeExtremeSingularValues()
122: @*/
123: PetscErrorCode  KSPMonitorSingularValue(KSP ksp,PetscInt n,PetscReal rnorm,PetscViewerAndFormat *dummy)
124: {
125:   PetscReal      emin,emax,c;
127:   PetscViewer    viewer = dummy->viewer;

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

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

150:    Collective on KSP

152:    Input Parameters:
153: +  ksp - the KSP context
154: .  its - iteration number
155: .  fgnorm - 2-norm of residual (or gradient)
156: -  dummy - a viewer

158:    Level: intermediate

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

164: .keywords: KSP, nonlinear, vector, monitor, view

166: .seealso: KSPMonitorSet(), KSPMonitorDefault(), VecView()
167: @*/
168: PetscErrorCode  KSPMonitorSolution(KSP ksp,PetscInt its,PetscReal fgnorm,PetscViewerAndFormat *dummy)
169: {
171:   Vec            x;
172:   PetscViewer    viewer = dummy->viewer;

176:   KSPBuildSolution(ksp,NULL,&x);
177:   PetscViewerPushFormat(viewer,dummy->format);
178:   VecView(x,viewer);
179:   PetscViewerPopFormat(viewer);
180:   return(0);
181: }

183: /*@C
184:    KSPMonitorDefault - Print the residual norm at each iteration of an
185:    iterative solver.

187:    Collective on KSP

189:    Input Parameters:
190: +  ksp   - iterative context
191: .  n     - iteration number
192: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).
193: -  dummy - an ASCII PetscViewer

195:    Level: intermediate

197: .keywords: KSP, default, monitor, residual

199: .seealso: KSPMonitorSet(), KSPMonitorTrueResidualNorm(), KSPMonitorLGResidualNormCreate()
200: @*/
201: PetscErrorCode  KSPMonitorDefault(KSP ksp,PetscInt n,PetscReal rnorm,PetscViewerAndFormat *dummy)
202: {
204:   PetscViewer    viewer =  dummy->viewer;

208:   PetscViewerPushFormat(viewer,dummy->format);
209:   PetscViewerASCIIAddTab(viewer,((PetscObject)ksp)->tablevel);
210:   if (n == 0 && ((PetscObject)ksp)->prefix) {
211:     PetscViewerASCIIPrintf(viewer,"  Residual norms for %s solve.\n",((PetscObject)ksp)->prefix);
212:   }
213:   PetscViewerASCIIPrintf(viewer,"%3D KSP Residual norm %14.12e \n",n,(double)rnorm);
214:   PetscViewerASCIISubtractTab(viewer,((PetscObject)ksp)->tablevel);
215:   PetscViewerPopFormat(viewer);
216:   return(0);
217: }

219: /*@C
220:    KSPMonitorTrueResidualNorm - Prints the true residual norm as well as the preconditioned
221:    residual norm at each iteration of an iterative solver.

223:    Collective on KSP

225:    Input Parameters:
226: +  ksp   - iterative context
227: .  n     - iteration number
228: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).
229: -  dummy - an ASCII PetscViewer

231:    Options Database Key:
232: .  -ksp_monitor_true_residual - Activates KSPMonitorTrueResidualNorm()

234:    Notes:
235:    When using right preconditioning, these values are equivalent.

237:    Level: intermediate

239: .keywords: KSP, default, monitor, residual

241: .seealso: KSPMonitorSet(), KSPMonitorDefault(), KSPMonitorLGResidualNormCreate(),KSPMonitorTrueResidualMaxNorm()
242: @*/
243: PetscErrorCode  KSPMonitorTrueResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,PetscViewerAndFormat *dummy)
244: {
246:   Vec            resid;
247:   PetscReal      truenorm,bnorm;
248:   PetscViewer    viewer = dummy->viewer;
249:   char           normtype[256];

253:   PetscViewerPushFormat(viewer,dummy->format);
254:   PetscViewerASCIIAddTab(viewer,((PetscObject)ksp)->tablevel);
255:   if (n == 0 && ((PetscObject)ksp)->prefix) {
256:     PetscViewerASCIIPrintf(viewer,"  Residual norms for %s solve.\n",((PetscObject)ksp)->prefix);
257:   }
258:   KSPBuildResidual(ksp,NULL,NULL,&resid);
259:   VecNorm(resid,NORM_2,&truenorm);
260:   VecDestroy(&resid);
261:   VecNorm(ksp->vec_rhs,NORM_2,&bnorm);
262:   PetscStrncpy(normtype,KSPNormTypes[ksp->normtype],sizeof(normtype));
263:   PetscStrtolower(normtype);
264:   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));
265:   PetscViewerASCIISubtractTab(viewer,((PetscObject)ksp)->tablevel);
266:   PetscViewerPopFormat(viewer);
267:   return(0);
268: }

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

273:    Collective on KSP

275:    Input Parameters:
276: +  ksp   - iterative context
277: .  n     - iteration number
278: .  rnorm - norm (preconditioned) residual value (may be estimated).
279: -  dummy - an ASCII viewer

281:    Options Database Key:
282: .  -ksp_monitor_max - Activates KSPMonitorTrueResidualMaxNorm()

284:    Notes:
285:    This could be implemented (better) with a flag in ksp.

287:    Level: intermediate

289: .keywords: KSP, default, monitor, residual

291: .seealso: KSPMonitorSet(), KSPMonitorDefault(), KSPMonitorLGResidualNormCreate(),KSPMonitorTrueResidualNorm()
292: @*/
293: PetscErrorCode  KSPMonitorTrueResidualMaxNorm(KSP ksp,PetscInt n,PetscReal rnorm,PetscViewerAndFormat *dummy)
294: {
296:   Vec            resid;
297:   PetscReal      truenorm,bnorm;
298:   PetscViewer    viewer = dummy->viewer;
299:   char           normtype[256];

303:   PetscViewerPushFormat(viewer,dummy->format);
304:   PetscViewerASCIIAddTab(viewer,((PetscObject)ksp)->tablevel);
305:   if (n == 0 && ((PetscObject)ksp)->prefix) {
306:     PetscViewerASCIIPrintf(viewer,"  Residual norms (max) for %s solve.\n",((PetscObject)ksp)->prefix);
307:   }
308:   KSPBuildResidual(ksp,NULL,NULL,&resid);
309:   VecNorm(resid,NORM_INFINITY,&truenorm);
310:   VecDestroy(&resid);
311:   VecNorm(ksp->vec_rhs,NORM_INFINITY,&bnorm);
312:   PetscStrncpy(normtype,KSPNormTypes[ksp->normtype],sizeof(normtype));
313:   PetscStrtolower(normtype);
314:   PetscViewerASCIIPrintf(viewer,"%3D KSP true resid max norm %14.12e ||r(i)||/||b|| %14.12e\n",n,(double)truenorm,(double)(truenorm/bnorm));
315:   PetscViewerASCIISubtractTab(viewer,((PetscObject)ksp)->tablevel);
316:   PetscViewerPopFormat(viewer);
317:   return(0);
318: }

320: PetscErrorCode  KSPMonitorRange_Private(KSP ksp,PetscInt it,PetscReal *per)
321: {
323:   Vec            resid;
324:   PetscReal      rmax,pwork;
325:   PetscInt       i,n,N;
326:   const PetscScalar *r;

329:   KSPBuildResidual(ksp,NULL,NULL,&resid);
330:   VecNorm(resid,NORM_INFINITY,&rmax);
331:   VecGetLocalSize(resid,&n);
332:   VecGetSize(resid,&N);
333:   VecGetArrayRead(resid,&r);
334:   pwork = 0.0;
335:   for (i=0; i<n; i++) pwork += (PetscAbsScalar(r[i]) > .20*rmax);
336:   VecRestoreArrayRead(resid,&r);
337:   VecDestroy(&resid);
338:   MPIU_Allreduce(&pwork,per,1,MPIU_REAL,MPIU_SUM,PetscObjectComm((PetscObject)ksp));
339:   *per = *per/N;
340:   return(0);
341: }

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

346:    Collective on KSP

348:    Input Parameters:
349: +  ksp   - iterative context
350: .  it    - iteration number
351: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).
352: -  dummy - an ASCII viewer

354:    Options Database Key:
355: .  -ksp_monitor_range - Activates KSPMonitorRange()

357:    Level: intermediate

359: .keywords: KSP, default, monitor, residual

361: .seealso: KSPMonitorSet(), KSPMonitorDefault(), KSPMonitorLGResidualNormCreate()
362: @*/
363: PetscErrorCode  KSPMonitorRange(KSP ksp,PetscInt it,PetscReal rnorm,PetscViewerAndFormat *dummy)
364: {
365:   PetscErrorCode   ierr;
366:   PetscReal        perc,rel;
367:   PetscViewer      viewer = dummy->viewer;
368:   /* should be in a MonitorRangeContext */
369:   static PetscReal prev;

373:   PetscViewerPushFormat(viewer,dummy->format);
374:   PetscViewerASCIIAddTab(viewer,((PetscObject)ksp)->tablevel);
375:   if (!it) prev = rnorm;
376:   if (it == 0 && ((PetscObject)ksp)->prefix) {
377:     PetscViewerASCIIPrintf(viewer,"  Residual norms for %s solve.\n",((PetscObject)ksp)->prefix);
378:   }
379:   KSPMonitorRange_Private(ksp,it,&perc);

381:   rel  = (prev - rnorm)/prev;
382:   prev = rnorm;
383:   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));
384:   PetscViewerASCIISubtractTab(viewer,((PetscObject)ksp)->tablevel);
385:   PetscViewerPopFormat(viewer);
386:   return(0);
387: }

389: /*@C
390:    KSPMonitorDynamicTolerance - Recompute the inner tolerance in every
391:    outer iteration in an adaptive way.

393:    Collective on KSP

395:    Input Parameters:
396: +  ksp   - iterative context
397: .  n     - iteration number (not used)
398: .  fnorm - the current residual norm
399: .  dummy - some context as a C struct. fields:
400:              coef: a scaling coefficient. default 1.0. can be passed through
401:                    -sub_ksp_dynamic_tolerance_param
402:              bnrm: norm of the right-hand side. store it to avoid repeated calculation

404:    Notes:
405:    This may be useful for a flexibly preconditioner Krylov method to
406:    control the accuracy of the inner solves needed to gaurantee the
407:    convergence of the outer iterations.

409:    Level: advanced

411: .keywords: KSP, inner tolerance

413: .seealso: KSPMonitorDynamicToleranceDestroy()
414: @*/
415: PetscErrorCode KSPMonitorDynamicTolerance(KSP ksp,PetscInt its,PetscReal fnorm,void *dummy)
416: {
418:   PC             pc;
419:   PetscReal      outer_rtol, outer_abstol, outer_dtol, inner_rtol;
420:   PetscInt       outer_maxits,nksp,first,i;
421:   KSPDynTolCtx   *scale   = (KSPDynTolCtx*)dummy;
422:   KSP            kspinner = NULL, *subksp = NULL;

425:   KSPGetPC(ksp, &pc);

427:   /* compute inner_rtol */
428:   if (scale->bnrm < 0.0) {
429:     Vec b;
430:     KSPGetRhs(ksp, &b);
431:     VecNorm(b, NORM_2, &(scale->bnrm));
432:   }
433:   KSPGetTolerances(ksp, &outer_rtol, &outer_abstol, &outer_dtol, &outer_maxits);
434:   inner_rtol = PetscMin(scale->coef * scale->bnrm * outer_rtol / fnorm, 0.999);
435:   /*PetscPrintf(PETSC_COMM_WORLD, "        Inner rtol = %g\n", (double)inner_rtol);*/

437:   /* if pc is ksp */
438:   PCKSPGetKSP(pc, &kspinner);
439:   if (kspinner) {
440:     KSPSetTolerances(kspinner, inner_rtol, outer_abstol, outer_dtol, outer_maxits);
441:     return(0);
442:   }

444:   /* if pc is bjacobi */
445:   PCBJacobiGetSubKSP(pc, &nksp, &first, &subksp);
446:   if (subksp) {
447:     for (i=0; i<nksp; i++) {
448:       KSPSetTolerances(subksp[i], inner_rtol, outer_abstol, outer_dtol, outer_maxits);
449:     }
450:     return(0);
451:   }

453:   /* todo: dynamic tolerance may apply to other types of pc too */
454:   return(0);
455: }

457: /*
458:   Destroy the dummy context used in KSPMonitorDynamicTolerance()
459: */
460: PetscErrorCode KSPMonitorDynamicToleranceDestroy(void **dummy)
461: {

465:   PetscFree(*dummy);
466:   return(0);
467: }

469: /*
470:   Default (short) KSP Monitor, same as KSPMonitorDefault() except
471:   it prints fewer digits of the residual as the residual gets smaller.
472:   This is because the later digits are meaningless and are often
473:   different on different machines; by using this routine different
474:   machines will usually generate the same output.
475: */
476: PetscErrorCode  KSPMonitorDefaultShort(KSP ksp,PetscInt its,PetscReal fnorm,PetscViewerAndFormat *dummy)
477: {
479:   PetscViewer    viewer = dummy->viewer;

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

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

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

505:    Collective on KSP

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

513:    Returns:
514: .  reason - KSP_CONVERGED_ITERATING, KSP_CONVERGED_ITS

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

523:    Level: advanced

525: .keywords: KSP, default, convergence, residual

527: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPSetNormType()
528: @*/
529: PetscErrorCode  KSPConvergedSkip(KSP ksp,PetscInt n,PetscReal rnorm,KSPConvergedReason *reason,void *dummy)
530: {
534:   *reason = KSP_CONVERGED_ITERATING;
535:   if (n >= ksp->max_it) *reason = KSP_CONVERGED_ITS;
536:   return(0);
537: }


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

543:    Collective on KSP

545:    Output Parameter:
546: .  ctx - convergence context

548:    Level: intermediate

550: .keywords: KSP, default, convergence, residual

552: .seealso: KSPConvergedDefault(), KSPConvergedDefaultDestroy(), KSPSetConvergenceTest(), KSPSetTolerances(),
553:           KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUIRNorm(), KSPConvergedDefaultSetUMIRNorm()
554: @*/
555: PetscErrorCode  KSPConvergedDefaultCreate(void **ctx)
556: {
557:   PetscErrorCode         ierr;
558:   KSPConvergedDefaultCtx *cctx;

561:   PetscNew(&cctx);
562:   *ctx = cctx;
563:   return(0);
564: }

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

571:    Collective on KSP

573:    Input Parameters:
574: .  ksp   - iterative context

576:    Options Database:
577: .   -ksp_converged_use_initial_residual_norm

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

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

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

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


590:    Level: intermediate

592: .keywords: KSP, default, convergence, residual

594: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUMIRNorm()
595: @*/
596: PetscErrorCode  KSPConvergedDefaultSetUIRNorm(KSP ksp)
597: {
598:   KSPConvergedDefaultCtx *ctx = (KSPConvergedDefaultCtx*) ksp->cnvP;

602:   if (ksp->converged != KSPConvergedDefault) return(0);
603:   if (ctx->mininitialrtol) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONGSTATE,"Cannot use KSPConvergedDefaultSetUIRNorm() and KSPConvergedDefaultSetUMIRNorm() together");
604:   ctx->initialrtol = PETSC_TRUE;
605:   return(0);
606: }

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

613:    Collective on KSP

615:    Input Parameters:
616: .  ksp   - iterative context

618:    Options Database:
619: .   -ksp_converged_use_min_initial_residual_norm

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

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

626:    Level: intermediate

628: .keywords: KSP, default, convergence, residual

630: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUIRNorm()
631: @*/
632: PetscErrorCode  KSPConvergedDefaultSetUMIRNorm(KSP ksp)
633: {
634:   KSPConvergedDefaultCtx *ctx = (KSPConvergedDefaultCtx*) ksp->cnvP;

638:   if (ksp->converged != KSPConvergedDefault) return(0);
639:   if (ctx->initialrtol) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONGSTATE,"Cannot use KSPConvergedDefaultSetUIRNorm() and KSPConvergedDefaultSetUMIRNorm() together");
640:   ctx->mininitialrtol = PETSC_TRUE;
641:   return(0);
642: }

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

647:    Collective on KSP

649:    Input Parameters:
650: +  ksp   - iterative context
651: .  n     - iteration number
652: .  rnorm - residual norm (may be estimated, depending on the method may be the preconditioned residual norm)
653: -  ctx - convergence context which must be created by KSPConvergedDefaultCreate()

655:    Output Parameter:
656: +   positive - if the iteration has converged;
657: .   negative - if residual norm exceeds divergence threshold;
658: -   0 - otherwise.

660:    Notes:
661:    KSPConvergedDefault() reaches convergence when   rnorm < MAX (rtol * rnorm_0, abstol);
662:    Divergence is detected if  rnorm > dtol * rnorm_0,

664:    where:
665: +     rtol = relative tolerance,
666: .     abstol = absolute tolerance.
667: .     dtol = divergence tolerance,
668: -     rnorm_0 is the two norm of the right hand side (or the preconditioned norm, depending on what was set with
669:           KSPSetNormType(). When initial guess is non-zero you
670:           can call KSPConvergedDefaultSetUIRNorm() to use the norm of (b - A*(initial guess))
671:           as the starting point for relative norm convergence testing, that is as rnorm_0

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

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

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

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

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

683:    Level: intermediate

685: .keywords: KSP, default, convergence, residual

687: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(), KSPConvergedReason, KSPGetConvergedReason(),
688:           KSPConvergedDefaultSetUIRNorm(), KSPConvergedDefaultSetUMIRNorm(), KSPConvergedDefaultCreate(), KSPConvergedDefaultDestroy()
689: @*/
690: PetscErrorCode  KSPConvergedDefault(KSP ksp,PetscInt n,PetscReal rnorm,KSPConvergedReason *reason,void *ctx)
691: {
692:   PetscErrorCode         ierr;
693:   KSPConvergedDefaultCtx *cctx = (KSPConvergedDefaultCtx*) ctx;
694:   KSPNormType            normtype;

699:   *reason = KSP_CONVERGED_ITERATING;

701:   KSPGetNormType(ksp,&normtype);
702:   if (normtype == KSP_NORM_NONE) return(0);

704:   if (!cctx) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_NULL,"Convergence context must have been created with KSPConvergedDefaultCreate()");
705:   if (!n) {
706:     /* if user gives initial guess need to compute norm of b */
707:     if (!ksp->guess_zero && !cctx->initialrtol) {
708:       PetscReal snorm = 0.0;
709:       if (ksp->normtype == KSP_NORM_UNPRECONDITIONED || ksp->pc_side == PC_RIGHT) {
710:         PetscInfo(ksp,"user has provided nonzero initial guess, computing 2-norm of RHS\n");
711:         VecNorm(ksp->vec_rhs,NORM_2,&snorm);        /*     <- b'*b */
712:       } else {
713:         Vec z;
714:         /* Should avoid allocating the z vector each time but cannot stash it in cctx because if KSPReset() is called the vector size might change */
715:         VecDuplicate(ksp->vec_rhs,&z);
716:         KSP_PCApply(ksp,ksp->vec_rhs,z);
717:         if (ksp->normtype == KSP_NORM_PRECONDITIONED) {
718:           PetscInfo(ksp,"user has provided nonzero initial guess, computing 2-norm of preconditioned RHS\n");
719:           VecNorm(z,NORM_2,&snorm);                 /*    dp <- b'*B'*B*b */
720:         } else if (ksp->normtype == KSP_NORM_NATURAL) {
721:           PetscScalar norm;
722:           PetscInfo(ksp,"user has provided nonzero initial guess, computing natural norm of RHS\n");
723:           VecDot(ksp->vec_rhs,z,&norm);
724:           snorm = PetscSqrtReal(PetscAbsScalar(norm));                            /*    dp <- b'*B*b */
725:         }
726:         VecDestroy(&z);
727:       }
728:       /* handle special case of zero RHS and nonzero guess */
729:       if (!snorm) {
730:         PetscInfo(ksp,"Special case, user has provided nonzero initial guess and zero RHS\n");
731:         snorm = rnorm;
732:       }
733:       if (cctx->mininitialrtol) ksp->rnorm0 = PetscMin(snorm,rnorm);
734:       else ksp->rnorm0 = snorm;
735:     } else {
736:       ksp->rnorm0 = rnorm;
737:     }
738:     ksp->ttol = PetscMax(ksp->rtol*ksp->rnorm0,ksp->abstol);
739:   }

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

743:   if (PetscIsInfOrNanReal(rnorm)) {
744:     PCFailedReason pcreason;
745:     PetscInt       sendbuf,pcreason_max;
746:     PCGetSetUpFailedReason(ksp->pc,&pcreason);
747:     sendbuf = (PetscInt)pcreason;
748:     MPI_Allreduce(&sendbuf,&pcreason_max,1,MPIU_INT,MPIU_MAX,PetscObjectComm((PetscObject)ksp));
749:     if (pcreason_max) {
750:       *reason = KSP_DIVERGED_PCSETUP_FAILED;
751:       VecSetInf(ksp->vec_sol);
752:       PetscInfo(ksp,"Linear solver pcsetup fails, declaring divergence \n");
753:     } else {
754:       *reason = KSP_DIVERGED_NANORINF;
755:       PetscInfo(ksp,"Linear solver has created a not a number (NaN) as the residual norm, declaring divergence \n");
756:     }
757:   } else if (rnorm <= ksp->ttol) {
758:     if (rnorm < ksp->abstol) {
759:       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);
760:       *reason = KSP_CONVERGED_ATOL;
761:     } else {
762:       if (cctx->initialrtol) {
763:         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);
764:       } else {
765:         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);
766:       }
767:       *reason = KSP_CONVERGED_RTOL;
768:     }
769:   } else if (rnorm >= ksp->divtol*ksp->rnorm0) {
770:     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);
771:     *reason = KSP_DIVERGED_DTOL;
772:   }
773:   return(0);
774: }

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

779:    Collective on KSP

781:    Input Parameters:
782: .  ctx - convergence context

784:    Level: intermediate

786: .keywords: KSP, default, convergence, residual

788: .seealso: KSPConvergedDefault(), KSPConvergedDefaultCreate(), KSPSetConvergenceTest(), KSPSetTolerances(), KSPConvergedSkip(),
789:           KSPConvergedReason, KSPGetConvergedReason(), KSPConvergedDefaultSetUIRNorm(), KSPConvergedDefaultSetUMIRNorm()
790: @*/
791: PetscErrorCode  KSPConvergedDefaultDestroy(void *ctx)
792: {
793:   PetscErrorCode         ierr;
794:   KSPConvergedDefaultCtx *cctx = (KSPConvergedDefaultCtx*) ctx;

797:   VecDestroy(&cctx->work);
798:   PetscFree(ctx);
799:   return(0);
800: }

802: /*
803:    KSPBuildSolutionDefault - Default code to create/move the solution.

805:    Input Parameters:
806: +  ksp - iterative context
807: -  v   - pointer to the user's vector

809:    Output Parameter:
810: .  V - pointer to a vector containing the solution

812:    Level: advanced

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

816: .keywords:  KSP, build, solution, default

818: .seealso: KSPGetSolution(), KSPBuildResidualDefault()
819: */
820: PetscErrorCode KSPBuildSolutionDefault(KSP ksp,Vec v,Vec *V)
821: {

825:   if (ksp->pc_side == PC_RIGHT) {
826:     if (ksp->pc) {
827:       if (v) {
828:         KSP_PCApply(ksp,ksp->vec_sol,v); *V = v;
829:       } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Not working with right preconditioner");
830:     } else {
831:       if (v) {
832:         VecCopy(ksp->vec_sol,v); *V = v;
833:       } else *V = ksp->vec_sol;
834:     }
835:   } else if (ksp->pc_side == PC_SYMMETRIC) {
836:     if (ksp->pc) {
837:       if (ksp->transpose_solve) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Not working with symmetric preconditioner and transpose solve");
838:       if (v) {
839:         PCApplySymmetricRight(ksp->pc,ksp->vec_sol,v);
840:         *V = v;
841:       } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Not working with symmetric preconditioner");
842:     } else {
843:       if (v) {
844:         VecCopy(ksp->vec_sol,v); *V = v;
845:       } else *V = ksp->vec_sol;
846:     }
847:   } else {
848:     if (v) {
849:       VecCopy(ksp->vec_sol,v); *V = v;
850:     } else *V = ksp->vec_sol;
851:   }
852:   return(0);
853: }

855: /*
856:    KSPBuildResidualDefault - Default code to compute the residual.

858:    Input Parameters:
859: .  ksp - iterative context
860: .  t   - pointer to temporary vector
861: .  v   - pointer to user vector

863:    Output Parameter:
864: .  V - pointer to a vector containing the residual

866:    Level: advanced

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

870: .keywords:  KSP, build, residual, default

872: .seealso: KSPBuildSolutionDefault()
873: */
874: PetscErrorCode KSPBuildResidualDefault(KSP ksp,Vec t,Vec v,Vec *V)
875: {
877:   Mat            Amat,Pmat;

880:   if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
881:   PCGetOperators(ksp->pc,&Amat,&Pmat);
882:   KSPBuildSolution(ksp,t,NULL);
883:   KSP_MatMult(ksp,Amat,t,v);
884:   VecAYPX(v,-1.0,ksp->vec_rhs);
885:   *V   = v;
886:   return(0);
887: }

889: /*@C
890:   KSPCreateVecs - Gets a number of work vectors.

892:   Input Parameters:
893: + ksp  - iterative context
894: . rightn  - number of right work vectors
895: - leftn   - number of left work vectors to allocate

897:   Output Parameter:
898: +  right - the array of vectors created
899: -  left - the array of left vectors

901:    Note: The right vector has as many elements as the matrix has columns. The left
902:      vector has as many elements as the matrix has rows.

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

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

909:    Level: advanced

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

913: @*/
914: PetscErrorCode KSPCreateVecs(KSP ksp,PetscInt rightn, Vec **right,PetscInt leftn,Vec **left)
915: {
917:   Vec            vecr = NULL,vecl = NULL;
918:   PetscBool      matset,pmatset;
919:   Mat            mat = NULL;

922:   if (rightn) {
923:     if (!right) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"You asked for right vectors but did not pass a pointer to hold them");
924:     if (ksp->vec_sol) vecr = ksp->vec_sol;
925:     else {
926:       if (ksp->pc) {
927:         PCGetOperatorsSet(ksp->pc,&matset,&pmatset);
928:         /* check for mat before pmat because for KSPLSQR pmat may be a different size than mat since pmat maybe mat'*mat */
929:         if (matset) {
930:           PCGetOperators(ksp->pc,&mat,NULL);
931:           MatCreateVecs(mat,&vecr,NULL);
932:         } else if (pmatset) {
933:           PCGetOperators(ksp->pc,NULL,&mat);
934:           MatCreateVecs(mat,&vecr,NULL);
935:         }
936:       }
937:       if (!vecr) {
938:         if (ksp->dm) {
939:           DMGetGlobalVector(ksp->dm,&vecr);
940:         } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONGSTATE,"You requested a vector from a KSP that cannot provide one");
941:       }
942:     }
943:     VecDuplicateVecs(vecr,rightn,right);
944:     if (!ksp->vec_sol) {
945:       if (mat) {
946:         VecDestroy(&vecr);
947:       } else if (ksp->dm) {
948:         DMRestoreGlobalVector(ksp->dm,&vecr);
949:       }
950:     }
951:   }
952:   if (leftn) {
953:     if (!left) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"You asked for left vectors but did not pass a pointer to hold them");
954:     if (ksp->vec_rhs) vecl = ksp->vec_rhs;
955:     else {
956:       if (ksp->pc) {
957:         PCGetOperatorsSet(ksp->pc,&matset,&pmatset);
958:         /* check for mat before pmat because for KSPLSQR pmat may be a different size than mat since pmat maybe mat'*mat */
959:         if (matset) {
960:           PCGetOperators(ksp->pc,&mat,NULL);
961:           MatCreateVecs(mat,NULL,&vecl);
962:         } else if (pmatset) {
963:           PCGetOperators(ksp->pc,NULL,&mat);
964:           MatCreateVecs(mat,NULL,&vecl);
965:         }
966:       }
967:       if (!vecl) {
968:         if (ksp->dm) {
969:           DMGetGlobalVector(ksp->dm,&vecl);
970:         } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONGSTATE,"You requested a vector from a KSP that cannot provide one");
971:       }
972:     }
973:     VecDuplicateVecs(vecl,leftn,left);
974:     if (!ksp->vec_rhs) {
975:       if (mat) {
976:         VecDestroy(&vecl);
977:       } else if (ksp->dm) {
978:         DMRestoreGlobalVector(ksp->dm,&vecl);
979:       }
980:     }
981:   }
982:   return(0);
983: }

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

988:   Collective on KSP

990:   Input Parameters:
991: + ksp  - iterative context
992: - nw   - number of work vectors to allocate

994:   Level: developer

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

998: @*/
999: PetscErrorCode KSPSetWorkVecs(KSP ksp,PetscInt nw)
1000: {

1004:   VecDestroyVecs(ksp->nwork,&ksp->work);
1005:   ksp->nwork = nw;
1006:   KSPCreateVecs(ksp,nw,&ksp->work,0,NULL);
1007:   PetscLogObjectParents(ksp,nw,ksp->work);
1008:   return(0);
1009: }

1011: /*
1012:   KSPDestroyDefault - Destroys a iterative context variable for methods with
1013:   no separate context.  Preferred calling sequence KSPDestroy().

1015:   Input Parameter:
1016: . ksp - the iterative context

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

1020: */
1021: PetscErrorCode KSPDestroyDefault(KSP ksp)
1022: {

1027:   PetscFree(ksp->data);
1028:   return(0);
1029: }

1031: /*@
1032:    KSPGetConvergedReason - Gets the reason the KSP iteration was stopped.

1034:    Not Collective

1036:    Input Parameter:
1037: .  ksp - the KSP context

1039:    Output Parameter:
1040: .  reason - negative value indicates diverged, positive value converged, see KSPConvergedReason

1042:    Possible values for reason: See also manual page for each reason
1043: $  KSP_CONVERGED_RTOL (residual 2-norm decreased by a factor of rtol, from 2-norm of right hand side)
1044: $  KSP_CONVERGED_ATOL (residual 2-norm less than abstol)
1045: $  KSP_CONVERGED_ITS (used by the preonly preconditioner that always uses ONE iteration, or when the KSPConvergedSkip() convergence test routine is set.
1046: $  KSP_CONVERGED_CG_NEG_CURVE (see note below)
1047: $  KSP_CONVERGED_CG_CONSTRAINED (see note below)
1048: $  KSP_CONVERGED_STEP_LENGTH (see note below)
1049: $  KSP_CONVERGED_ITERATING (returned if the solver is not yet finished)
1050: $  KSP_DIVERGED_ITS  (required more than its to reach convergence)
1051: $  KSP_DIVERGED_DTOL (residual norm increased by a factor of divtol)
1052: $  KSP_DIVERGED_NANORINF (residual norm became Not-a-number or Inf likely due to 0/0)
1053: $  KSP_DIVERGED_BREAKDOWN (generic breakdown in method)
1054: $  KSP_DIVERGED_BREAKDOWN_BICG (Initial residual is orthogonal to preconditioned initial residual. Try a different preconditioner, or a different initial Level.)

1056:    Options Database:
1057: .   -ksp_converged_reason - prints the reason to standard out

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

1061:    The values  KSP_CONVERGED_CG_NEG_CURVE, KSP_CONVERGED_CG_CONSTRAINED, and KSP_CONVERGED_STEP_LENGTH are returned only by the special KSPCGNASH, KSPCGSTCG, and KSPCGGLTR
1062:    solvers which are used by the SNESNEWTONTR (trust region) solver.

1064:    Level: intermediate

1066: .keywords: KSP, nonlinear, set, convergence, test

1068: .seealso: KSPSetConvergenceTest(), KSPConvergedDefault(), KSPSetTolerances(), KSPConvergedReason
1069: @*/
1070: PetscErrorCode  KSPGetConvergedReason(KSP ksp,KSPConvergedReason *reason)
1071: {
1075:   *reason = ksp->reason;
1076:   return(0);
1077: }

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

1083:    Logically Collective on KSP

1085:    Input Parameters:
1086: +  ksp - the preconditioner context
1087: -  dm - the dm, cannot be NULL

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

1093:    Level: intermediate

1095: .seealso: KSPGetDM(), KSPSetDMActive(), KSPSetComputeOperators(), KSPSetComputeRHS(), KSPSetComputeInitialGuess(), DMKSPSetComputeOperators(), DMKSPSetComputeRHS(), DMKSPSetComputeInitialGuess()
1096: @*/
1097: PetscErrorCode  KSPSetDM(KSP ksp,DM dm)
1098: {
1100:   PC             pc;

1105:   PetscObjectReference((PetscObject)dm);
1106:   if (ksp->dm) {                /* Move the DMSNES context over to the new DM unless the new DM already has one */
1107:     if (ksp->dm->dmksp && !dm->dmksp) {
1108:       DMKSP kdm;
1109:       DMCopyDMKSP(ksp->dm,dm);
1110:       DMGetDMKSP(ksp->dm,&kdm);
1111:       if (kdm->originaldm == ksp->dm) kdm->originaldm = dm; /* Grant write privileges to the replacement DM */
1112:     }
1113:     DMDestroy(&ksp->dm);
1114:   }
1115:   ksp->dm       = dm;
1116:   ksp->dmAuto   = PETSC_FALSE;
1117:   KSPGetPC(ksp,&pc);
1118:   PCSetDM(pc,dm);
1119:   ksp->dmActive = PETSC_TRUE;
1120:   return(0);
1121: }

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

1126:    Logically Collective on KSP

1128:    Input Parameters:
1129: +  ksp - the preconditioner context
1130: -  flg - use the DM

1132:    Notes:
1133:    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.

1135:    Level: intermediate

1137: .seealso: KSPGetDM(), KSPSetDM(), SNESSetDM(), KSPSetComputeOperators(), KSPSetComputeRHS(), KSPSetComputeInitialGuess()
1138: @*/
1139: PetscErrorCode  KSPSetDMActive(KSP ksp,PetscBool flg)
1140: {
1144:   ksp->dmActive = flg;
1145:   return(0);
1146: }

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

1151:    Not Collective

1153:    Input Parameter:
1154: . ksp - the preconditioner context

1156:    Output Parameter:
1157: .  dm - the dm

1159:    Level: intermediate


1162: .seealso: KSPSetDM(), KSPSetDMActive()
1163: @*/
1164: PetscErrorCode  KSPGetDM(KSP ksp,DM *dm)
1165: {

1170:   if (!ksp->dm) {
1171:     DMShellCreate(PetscObjectComm((PetscObject)ksp),&ksp->dm);
1172:     ksp->dmAuto = PETSC_TRUE;
1173:   }
1174:   *dm = ksp->dm;
1175:   return(0);
1176: }

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

1181:    Logically Collective on KSP

1183:    Input Parameters:
1184: +  ksp - the KSP context
1185: -  usrP - optional user context

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

1190:    Level: intermediate

1192: .keywords: KSP, set, application, context

1194: .seealso: KSPGetApplicationContext()
1195: @*/
1196: PetscErrorCode  KSPSetApplicationContext(KSP ksp,void *usrP)
1197: {
1199:   PC             pc;

1203:   ksp->user = usrP;
1204:   KSPGetPC(ksp,&pc);
1205:   PCSetApplicationContext(pc,usrP);
1206:   return(0);
1207: }

1209: /*@
1210:    KSPGetApplicationContext - Gets the user-defined context for the linear solver.

1212:    Not Collective

1214:    Input Parameter:
1215: .  ksp - KSP context

1217:    Output Parameter:
1218: .  usrP - user context

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

1223:    Level: intermediate

1225: .keywords: KSP, get, application, context

1227: .seealso: KSPSetApplicationContext()
1228: @*/
1229: PetscErrorCode  KSPGetApplicationContext(KSP ksp,void *usrP)
1230: {
1233:   *(void**)usrP = ksp->user;
1234:   return(0);
1235: }