Actual source code: kspimpl.h


  2: #ifndef _KSPIMPL_H
  3:   #define _KSPIMPL_H

  5: #include <petscksp.h>
  6: #include <petscds.h>
  7: #include <petsc/private/petscimpl.h>

  9: /* SUBMANSEC = KSP */

 11: PETSC_EXTERN PetscBool      KSPRegisterAllCalled;
 12: PETSC_EXTERN PetscBool      KSPMonitorRegisterAllCalled;
 13: PETSC_EXTERN PetscErrorCode KSPRegisterAll(void);
 14: PETSC_EXTERN PetscErrorCode KSPMonitorRegisterAll(void);
 15: PETSC_EXTERN PetscErrorCode KSPGuessRegisterAll(void);
 16: PETSC_EXTERN PetscErrorCode KSPMatRegisterAll(void);

 18: typedef struct _KSPOps *KSPOps;

 20: struct _KSPOps {
 21:   PetscErrorCode (*buildsolution)(KSP, Vec, Vec *);      /* Returns a pointer to the solution, or
 22:                                                           calculates the solution in a
 23:                                                           user-provided area. */
 24:   PetscErrorCode (*buildresidual)(KSP, Vec, Vec, Vec *); /* Returns a pointer to the residual, or
 25:                                                           calculates the residual in a
 26:                                                           user-provided area.  */
 27:   PetscErrorCode (*solve)(KSP);                          /* actual solver */
 28:   PetscErrorCode (*matsolve)(KSP, Mat, Mat);             /* multiple dense RHS solver */
 29:   PetscErrorCode (*setup)(KSP);
 30:   PetscErrorCode (*setfromoptions)(KSP, PetscOptionItems *);
 31:   PetscErrorCode (*publishoptions)(KSP);
 32:   PetscErrorCode (*computeextremesingularvalues)(KSP, PetscReal *, PetscReal *);
 33:   PetscErrorCode (*computeeigenvalues)(KSP, PetscInt, PetscReal *, PetscReal *, PetscInt *);
 34:   PetscErrorCode (*computeritz)(KSP, PetscBool, PetscBool, PetscInt *, Vec[], PetscReal *, PetscReal *);
 35:   PetscErrorCode (*destroy)(KSP);
 36:   PetscErrorCode (*view)(KSP, PetscViewer);
 37:   PetscErrorCode (*reset)(KSP);
 38:   PetscErrorCode (*load)(KSP, PetscViewer);
 39: };

 41: typedef struct _KSPGuessOps *KSPGuessOps;

 43: struct _KSPGuessOps {
 44:   PetscErrorCode (*formguess)(KSPGuess, Vec, Vec); /* Form initial guess */
 45:   PetscErrorCode (*update)(KSPGuess, Vec, Vec);    /* Update database */
 46:   PetscErrorCode (*setfromoptions)(KSPGuess);
 47:   PetscErrorCode (*settolerance)(KSPGuess, PetscReal);
 48:   PetscErrorCode (*setup)(KSPGuess);
 49:   PetscErrorCode (*destroy)(KSPGuess);
 50:   PetscErrorCode (*view)(KSPGuess, PetscViewer);
 51:   PetscErrorCode (*reset)(KSPGuess);
 52: };

 54: /*
 55:    Defines the KSPGuess data structure.
 56: */
 57: struct _p_KSPGuess {
 58:   PETSCHEADER(struct _KSPGuessOps);
 59:   KSP              ksp;       /* the parent KSP */
 60:   Mat              A;         /* the current linear operator */
 61:   PetscObjectState omatstate; /* previous linear operator state */
 62:   void            *data;      /* pointer to the specific implementation */
 63: };

 65: PETSC_EXTERN PetscErrorCode KSPGuessCreate_Fischer(KSPGuess);
 66: PETSC_EXTERN PetscErrorCode KSPGuessCreate_POD(KSPGuess);

 68:   /*
 69:      Maximum number of monitors you can run with a single KSP
 70: */
 71:   #define MAXKSPMONITORS    5
 72:   #define MAXKSPREASONVIEWS 5
 73: typedef enum {
 74:   KSP_SETUP_NEW,
 75:   KSP_SETUP_NEWMATRIX,
 76:   KSP_SETUP_NEWRHS
 77: } KSPSetUpStage;

 79: /*
 80:    Defines the KSP data structure.
 81: */
 82: struct _p_KSP {
 83:   PETSCHEADER(struct _KSPOps);
 84:   DM        dm;
 85:   PetscBool dmAuto;   /* DM was created automatically by KSP */
 86:   PetscBool dmActive; /* KSP should use DM for computing operators */
 87:   /*------------------------- User parameters--------------------------*/
 88:   PetscInt  max_it; /* maximum number of iterations */
 89:   KSPGuess  guess;
 90:   PetscBool guess_zero,                                  /* flag for whether initial guess is 0 */
 91:     calc_sings,                                          /* calculate extreme Singular Values */
 92:     calc_ritz,                                           /* calculate (harmonic) Ritz pairs */
 93:     guess_knoll;                                         /* use initial guess of PCApply(ksp->B,b */
 94:   PCSide    pc_side;                                     /* flag for left, right, or symmetric preconditioning */
 95:   PetscInt  normsupporttable[KSP_NORM_MAX][PC_SIDE_MAX]; /* Table of supported norms and pc_side, see KSPSetSupportedNorm() */
 96:   PetscReal rtol,                                        /* relative tolerance */
 97:     abstol,                                              /* absolute tolerance */
 98:     ttol,                                                /* (not set by user)  */
 99:     divtol;                                              /* divergence tolerance */
100:   PetscReal          rnorm0;                             /* initial residual norm (used for divergence testing) */
101:   PetscReal          rnorm;                              /* current residual norm */
102:   KSPConvergedReason reason;
103:   PetscBool          errorifnotconverged; /* create an error if the KSPSolve() does not converge */

105:   Vec        vec_sol, vec_rhs; /* pointer to where user has stashed
106:                                       the solution and rhs, these are
107:                                       never touched by the code, only
108:                                       passed back to the user */
109:   PetscReal *res_hist;         /* If !0 stores residual each at iteration */
110:   PetscReal *res_hist_alloc;   /* If !0 means user did not provide buffer, needs deallocation */
111:   size_t     res_hist_len;     /* current size of residual history array */
112:   size_t     res_hist_max;     /* actual amount of storage in residual history */
113:   PetscBool  res_hist_reset;   /* reset history to length zero for each new solve */
114:   PetscReal *err_hist;         /* If !0 stores error at each iteration */
115:   PetscReal *err_hist_alloc;   /* If !0 means user did not provide buffer, needs deallocation */
116:   size_t     err_hist_len;     /* current size of error history array */
117:   size_t     err_hist_max;     /* actual amount of storage in error history */
118:   PetscBool  err_hist_reset;   /* reset history to length zero for each new solve */

120:   PetscInt  chknorm; /* only compute/check norm if iterations is great than this */
121:   PetscBool lagnorm; /* Lag the residual norm calculation so that it is computed as part of the
122:                                         MPI_Allreduce() for computing the inner products for the next iteration. */

124:   PetscInt nmax; /* maximum number of right-hand sides to be handled simultaneously */

126:   /* --------User (or default) routines (most return -1 on error) --------*/
127:   PetscErrorCode (*monitor[MAXKSPMONITORS])(KSP, PetscInt, PetscReal, void *); /* returns control to user after */
128:   PetscErrorCode (*monitordestroy[MAXKSPMONITORS])(void **);                   /* */
129:   void     *monitorcontext[MAXKSPMONITORS];                                    /* residual calculation, allows user */
130:   PetscInt  numbermonitors;                                                    /* to, for instance, print residual norm, etc. */
131:   PetscBool pauseFinal;                                                        /* Pause all drawing monitor at the final iterate */

133:   PetscErrorCode (*reasonview[MAXKSPREASONVIEWS])(KSP, void *);    /* KSP converged reason view */
134:   PetscErrorCode (*reasonviewdestroy[MAXKSPREASONVIEWS])(void **); /* Optional destroy routine */
135:   void    *reasonviewcontext[MAXKSPREASONVIEWS];                   /* User context */
136:   PetscInt numberreasonviews;                                      /* Number if reason viewers */

138:   PetscErrorCode (*converged)(KSP, PetscInt, PetscReal, KSPConvergedReason *, void *);
139:   PetscErrorCode (*convergeddestroy)(void *);
140:   void *cnvP;

142:   void *user; /* optional user-defined context */

144:   PC pc;

146:   void *data; /* holder for misc stuff associated
147:                                    with a particular iterative solver */

149:   PetscBool         view, viewPre, viewRate, viewMat, viewPMat, viewRhs, viewSol, viewMatExp, viewEV, viewSV, viewEVExp, viewFinalRes, viewPOpExp, viewDScale;
150:   PetscViewer       viewer, viewerPre, viewerRate, viewerMat, viewerPMat, viewerRhs, viewerSol, viewerMatExp, viewerEV, viewerSV, viewerEVExp, viewerFinalRes, viewerPOpExp, viewerDScale;
151:   PetscViewerFormat format, formatPre, formatRate, formatMat, formatPMat, formatRhs, formatSol, formatMatExp, formatEV, formatSV, formatEVExp, formatFinalRes, formatPOpExp, formatDScale;

153:   /* ----------------Default work-area management -------------------- */
154:   PetscInt nwork;
155:   Vec     *work;

157:   KSPSetUpStage setupstage;
158:   PetscBool     setupnewmatrix; /* true if we need to call ksp->ops->setup with KSP_SETUP_NEWMATRIX */

160:   PetscInt its;      /* number of iterations so far computed in THIS linear solve*/
161:   PetscInt totalits; /* number of iterations used by this KSP object since it was created */

163:   PetscBool transpose_solve; /* solve transpose system instead */
164:   struct {
165:     Mat       AT, BT;
166:     PetscBool use_explicittranspose; /* transpose the system explicitly in KSPSolveTranspose */
167:     PetscBool reuse_transpose;       /* reuse the previous transposed system */
168:   } transpose;

170:   KSPNormType normtype; /* type of norm used for convergence tests */

172:   PCSide      pc_side_set;  /* PC type set explicitly by user */
173:   KSPNormType normtype_set; /* Norm type set explicitly by user */

175:   /*   Allow diagonally scaling the matrix before computing the preconditioner or using
176:        the Krylov method. Note this is NOT just Jacobi preconditioning */

178:   PetscBool dscale;     /* diagonal scale system; used with KSPSetDiagonalScale() */
179:   PetscBool dscalefix;  /* unscale system after solve */
180:   PetscBool dscalefix2; /* system has been unscaled */
181:   Vec       diagonal;   /* 1/sqrt(diag of matrix) */
182:   Vec       truediagonal;

184:   /* Allow declaring convergence when negative curvature is detected */
185:   PetscBool converged_neg_curve;

187:   PetscInt  setfromoptionscalled;
188:   PetscBool skippcsetfromoptions; /* if set then KSPSetFromOptions() does not call PCSetFromOptions() */

190:   PetscErrorCode (*presolve)(KSP, Vec, Vec, void *);
191:   PetscErrorCode (*postsolve)(KSP, Vec, Vec, void *);
192:   void *prectx, *postctx;
193: };

195: typedef struct { /* dummy data structure used in KSPMonitorDynamicTolerance() */
196:   PetscReal coef;
197:   PetscReal bnrm;
198: } KSPDynTolCtx;

200: typedef struct {
201:   PetscBool initialrtol;    /* default relative residual decrease is computed from initial residual, not rhs */
202:   PetscBool mininitialrtol; /* default relative residual decrease is computed from min of initial residual and rhs */
203:   PetscBool convmaxits;     /* if true, the convergence test returns KSP_CONVERGED_ITS if the maximum number of iterations is reached */
204:   Vec       work;
205: } KSPConvergedDefaultCtx;

207: static inline PetscErrorCode KSPLogResidualHistory(KSP ksp, PetscReal norm)
208: {
209:   PetscFunctionBegin;
210:   PetscCall(PetscObjectSAWsTakeAccess((PetscObject)ksp));
211:   if (ksp->res_hist && ksp->res_hist_max > ksp->res_hist_len) ksp->res_hist[ksp->res_hist_len++] = norm;
212:   PetscCall(PetscObjectSAWsGrantAccess((PetscObject)ksp));
213:   PetscFunctionReturn(PETSC_SUCCESS);
214: }

216: static inline PetscErrorCode KSPLogErrorHistory(KSP ksp)
217: {
218:   DM dm;

220:   PetscFunctionBegin;
221:   PetscCall(PetscObjectSAWsTakeAccess((PetscObject)ksp));
222:   PetscCall(KSPGetDM(ksp, &dm));
223:   if (dm && ksp->err_hist && ksp->err_hist_max > ksp->err_hist_len) {
224:     PetscSimplePointFunc exactSol;
225:     void                *exactCtx;
226:     PetscDS              ds;
227:     Vec                  u;
228:     PetscReal            error;
229:     PetscInt             Nf;

231:     PetscCall(KSPBuildSolution(ksp, NULL, &u));
232:     /* TODO Was needed to correct for Newton solution, but I just need to set a solution */
233:     //PetscCall(VecScale(u, -1.0));
234:     /* TODO Case when I have a solution */
235:     if (0) {
236:       PetscCall(DMGetDS(dm, &ds));
237:       PetscCall(PetscDSGetNumFields(ds, &Nf));
238:       PetscCheck(Nf <= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cannot handle number of fields %" PetscInt_FMT " > 1 right now", Nf);
239:       PetscCall(PetscDSGetExactSolution(ds, 0, &exactSol, &exactCtx));
240:       PetscCall(DMComputeL2FieldDiff(dm, 0.0, &exactSol, &exactCtx, u, &error));
241:     } else {
242:       /* The null solution A 0 = 0 */
243:       PetscCall(VecNorm(u, NORM_2, &error));
244:     }
245:     ksp->err_hist[ksp->err_hist_len++] = error;
246:   }
247:   PetscCall(PetscObjectSAWsGrantAccess((PetscObject)ksp));
248:   PetscFunctionReturn(PETSC_SUCCESS);
249: }

251: static inline PetscScalar KSPNoisyHash_Private(PetscInt xx)
252: {
253:   unsigned int x = (unsigned int)xx;
254:   x              = ((x >> 16) ^ x) * 0x45d9f3b;
255:   x              = ((x >> 16) ^ x) * 0x45d9f3b;
256:   x              = ((x >> 16) ^ x);
257:   return (PetscScalar)((PetscInt64)x - 2147483648) * 5.e-10; /* center around zero, scaled about -1. to 1.*/
258: }

260: static inline PetscErrorCode KSPSetNoisy_Private(Vec v)
261: {
262:   PetscScalar *a;
263:   PetscInt     n, istart;

265:   PetscFunctionBegin;
266:   PetscCall(VecGetOwnershipRange(v, &istart, NULL));
267:   PetscCall(VecGetLocalSize(v, &n));
268:   PetscCall(VecGetArrayWrite(v, &a));
269:   for (PetscInt i = 0; i < n; ++i) a[i] = KSPNoisyHash_Private(i + istart);
270:   PetscCall(VecRestoreArrayWrite(v, &a));
271:   PetscFunctionReturn(PETSC_SUCCESS);
272: }

274: PETSC_INTERN PetscErrorCode KSPSetUpNorms_Private(KSP, PetscBool, KSPNormType *, PCSide *);

276: PETSC_INTERN PetscErrorCode KSPPlotEigenContours_Private(KSP, PetscInt, const PetscReal *, const PetscReal *);

278: typedef struct _p_DMKSP  *DMKSP;
279: typedef struct _DMKSPOps *DMKSPOps;
280: struct _DMKSPOps {
281:   PetscErrorCode (*computeoperators)(KSP, Mat, Mat, void *);
282:   PetscErrorCode (*computerhs)(KSP, Vec, void *);
283:   PetscErrorCode (*computeinitialguess)(KSP, Vec, void *);
284:   PetscErrorCode (*destroy)(DMKSP *);
285:   PetscErrorCode (*duplicate)(DMKSP, DMKSP);
286: };

288: struct _p_DMKSP {
289:   PETSCHEADER(struct _DMKSPOps);
290:   void *operatorsctx;
291:   void *rhsctx;
292:   void *initialguessctx;
293:   void *data;

295:   /* This is NOT reference counted. The DM on which this context was first created is cached here to implement one-way
296:    * copy-on-write. When DMGetDMKSPWrite() sees a request using a different DM, it makes a copy. Thus, if a user
297:    * only interacts directly with one level, e.g., using KSPSetComputeOperators(), then coarse levels are constructed by
298:    * PCMG, then the user changes the routine with another call to KSPSetComputeOperators(), it automatically propagates
299:    * to all the levels. If instead, they get out a specific level and set the routines on that level, subsequent changes
300:    * to the original level will no longer propagate to that level.
301:    */
302:   DM originaldm;

304:   void (*fortran_func_pointers[3])(void); /* Store our own function pointers so they are associated with the DMKSP instead of the DM */
305: };
306: PETSC_EXTERN PetscErrorCode DMGetDMKSP(DM, DMKSP *);
307: PETSC_EXTERN PetscErrorCode DMGetDMKSPWrite(DM, DMKSP *);
308: PETSC_EXTERN PetscErrorCode DMCopyDMKSP(DM, DM);

310: /*
311:        These allow the various Krylov methods to apply to either the linear system or its transpose.
312: */
313: static inline PetscErrorCode KSP_RemoveNullSpace(KSP ksp, Vec y)
314: {
315:   PetscFunctionBegin;
316:   if (ksp->pc_side == PC_LEFT) {
317:     Mat          A;
318:     MatNullSpace nullsp;

320:     PetscCall(PCGetOperators(ksp->pc, &A, NULL));
321:     PetscCall(MatGetNullSpace(A, &nullsp));
322:     if (nullsp) PetscCall(MatNullSpaceRemove(nullsp, y));
323:   }
324:   PetscFunctionReturn(PETSC_SUCCESS);
325: }

327: static inline PetscErrorCode KSP_RemoveNullSpaceTranspose(KSP ksp, Vec y)
328: {
329:   PetscFunctionBegin;
330:   if (ksp->pc_side == PC_LEFT) {
331:     Mat          A;
332:     MatNullSpace nullsp;

334:     PetscCall(PCGetOperators(ksp->pc, &A, NULL));
335:     PetscCall(MatGetTransposeNullSpace(A, &nullsp));
336:     if (nullsp) PetscCall(MatNullSpaceRemove(nullsp, y));
337:   }
338:   PetscFunctionReturn(PETSC_SUCCESS);
339: }

341: static inline PetscErrorCode KSP_MatMult(KSP ksp, Mat A, Vec x, Vec y)
342: {
343:   PetscFunctionBegin;
344:   if (ksp->transpose_solve) PetscCall(MatMultTranspose(A, x, y));
345:   else PetscCall(MatMult(A, x, y));
346:   PetscFunctionReturn(PETSC_SUCCESS);
347: }

349: static inline PetscErrorCode KSP_MatMultTranspose(KSP ksp, Mat A, Vec x, Vec y)
350: {
351:   PetscFunctionBegin;
352:   if (ksp->transpose_solve) PetscCall(MatMult(A, x, y));
353:   else PetscCall(MatMultTranspose(A, x, y));
354:   PetscFunctionReturn(PETSC_SUCCESS);
355: }

357: static inline PetscErrorCode KSP_MatMultHermitianTranspose(KSP ksp, Mat A, Vec x, Vec y)
358: {
359:   PetscFunctionBegin;
360:   if (!ksp->transpose_solve) PetscCall(MatMultHermitianTranspose(A, x, y));
361:   else {
362:     Vec w;

364:     PetscCall(VecDuplicate(x, &w));
365:     PetscCall(VecCopy(x, w));
366:     PetscCall(VecConjugate(w));
367:     PetscCall(MatMult(A, w, y));
368:     PetscCall(VecDestroy(&w));
369:     PetscCall(VecConjugate(y));
370:   }
371:   PetscFunctionReturn(PETSC_SUCCESS);
372: }

374: static inline PetscErrorCode KSP_PCApply(KSP ksp, Vec x, Vec y)
375: {
376:   PetscFunctionBegin;
377:   if (ksp->transpose_solve) {
378:     PetscCall(PCApplyTranspose(ksp->pc, x, y));
379:     PetscCall(KSP_RemoveNullSpaceTranspose(ksp, y));
380:   } else {
381:     PetscCall(PCApply(ksp->pc, x, y));
382:     PetscCall(KSP_RemoveNullSpace(ksp, y));
383:   }
384:   PetscFunctionReturn(PETSC_SUCCESS);
385: }

387: static inline PetscErrorCode KSP_PCApplyTranspose(KSP ksp, Vec x, Vec y)
388: {
389:   PetscFunctionBegin;
390:   if (ksp->transpose_solve) {
391:     PetscCall(PCApply(ksp->pc, x, y));
392:     PetscCall(KSP_RemoveNullSpace(ksp, y));
393:   } else {
394:     PetscCall(PCApplyTranspose(ksp->pc, x, y));
395:     PetscCall(KSP_RemoveNullSpaceTranspose(ksp, y));
396:   }
397:   PetscFunctionReturn(PETSC_SUCCESS);
398: }

400: static inline PetscErrorCode KSP_PCApplyHermitianTranspose(KSP ksp, Vec x, Vec y)
401: {
402:   PetscFunctionBegin;
403:   PetscCall(VecConjugate(x));
404:   PetscCall(KSP_PCApplyTranspose(ksp, x, y));
405:   PetscCall(VecConjugate(x));
406:   PetscCall(VecConjugate(y));
407:   PetscFunctionReturn(PETSC_SUCCESS);
408: }

410: static inline PetscErrorCode KSP_PCMatApply(KSP ksp, Mat X, Mat Y)
411: {
412:   PetscFunctionBegin;
413:   if (ksp->transpose_solve) {
414:     PetscBool flg;
415:     PetscCall(PetscObjectTypeCompareAny((PetscObject)ksp->pc, &flg, PCNONE, PCICC, PCCHOLESKY, ""));
416:     PetscCheck(flg, PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "PCMatApplyTranspose() not yet implemented for nonsymmetric PC");
417:   }
418:   PetscCall(PCMatApply(ksp->pc, X, Y));
419:   PetscFunctionReturn(PETSC_SUCCESS);
420: }

422: static inline PetscErrorCode KSP_PCMatApplyTranspose(KSP ksp, Mat X, Mat Y)
423: {
424:   PetscFunctionBegin;
425:   if (!ksp->transpose_solve) {
426:     PetscBool flg;
427:     PetscCall(PetscObjectTypeCompareAny((PetscObject)ksp->pc, &flg, PCNONE, PCICC, PCCHOLESKY, ""));
428:     PetscCheck(flg, PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "PCMatApplyTranspose() not yet implemented for nonsymmetric PC");
429:   }
430:   PetscCall(PCMatApply(ksp->pc, X, Y));
431:   PetscFunctionReturn(PETSC_SUCCESS);
432: }

434: static inline PetscErrorCode KSP_PCApplyBAorAB(KSP ksp, Vec x, Vec y, Vec w)
435: {
436:   PetscFunctionBegin;
437:   if (ksp->transpose_solve) {
438:     PetscCall(PCApplyBAorABTranspose(ksp->pc, ksp->pc_side, x, y, w));
439:     PetscCall(KSP_RemoveNullSpaceTranspose(ksp, y));
440:   } else {
441:     PetscCall(PCApplyBAorAB(ksp->pc, ksp->pc_side, x, y, w));
442:     PetscCall(KSP_RemoveNullSpace(ksp, y));
443:   }
444:   PetscFunctionReturn(PETSC_SUCCESS);
445: }

447: static inline PetscErrorCode KSP_PCApplyBAorABTranspose(KSP ksp, Vec x, Vec y, Vec w)
448: {
449:   PetscFunctionBegin;
450:   if (ksp->transpose_solve) PetscCall(PCApplyBAorAB(ksp->pc, ksp->pc_side, x, y, w));
451:   else PetscCall(PCApplyBAorABTranspose(ksp->pc, ksp->pc_side, x, y, w));
452:   PetscFunctionReturn(PETSC_SUCCESS);
453: }

455: PETSC_EXTERN PetscLogEvent KSP_GMRESOrthogonalization;
456: PETSC_EXTERN PetscLogEvent KSP_SetUp;
457: PETSC_EXTERN PetscLogEvent KSP_Solve;
458: PETSC_EXTERN PetscLogEvent KSP_Solve_FS_0;
459: PETSC_EXTERN PetscLogEvent KSP_Solve_FS_1;
460: PETSC_EXTERN PetscLogEvent KSP_Solve_FS_2;
461: PETSC_EXTERN PetscLogEvent KSP_Solve_FS_3;
462: PETSC_EXTERN PetscLogEvent KSP_Solve_FS_4;
463: PETSC_EXTERN PetscLogEvent KSP_Solve_FS_S;
464: PETSC_EXTERN PetscLogEvent KSP_Solve_FS_L;
465: PETSC_EXTERN PetscLogEvent KSP_Solve_FS_U;
466: PETSC_EXTERN PetscLogEvent KSP_SolveTranspose;
467: PETSC_EXTERN PetscLogEvent KSP_MatSolve;
468: PETSC_EXTERN PetscLogEvent KSP_MatSolveTranspose;

470: PETSC_INTERN PetscErrorCode MatGetSchurComplement_Basic(Mat, IS, IS, IS, IS, MatReuse, Mat *, MatSchurComplementAinvType, MatReuse, Mat *);
471: PETSC_INTERN PetscErrorCode PCPreSolveChangeRHS(PC, PetscBool *);

473:   /*MC
474:    KSPCheckDot - Checks if the result of a dot product used by the corresponding `KSP` contains Inf or NaN. These indicate that the previous
475:       application of the preconditioner generated an error. Sets a `KSPConvergedReason` and returns if the `PC` set a `PCFailedReason`.

477:    Collective

479:    Input Parameter:
480: .  ksp - the linear solver `KSP` context.

482:    Output Parameter:
483: .  beta - the result of the inner product

485:    Level: developer

487:    Developer Notes:
488:    Used to manage returning from `KSP` solvers whose preconditioners have failed, possibly only a subset of MPI ranks, in some way

490:    It uses the fact that `KSP` piggy-backs the collectivity of certain error conditions on the results of norms and inner products.

492: .seealso: `PCFailedReason`, `KSPConvergedReason`, `PCGetFailedReasonRank()`, `KSP`, `KSPCreate()`, `KSPSetType()`, `KSP`, `KSPCheckNorm()`, `KSPCheckSolve()`
493: M*/
494:   #define KSPCheckDot(ksp, beta) \
495:     do { \
496:       if (PetscIsInfOrNanScalar(beta)) { \
497:         PetscCheck(!ksp->errorifnotconverged, PetscObjectComm((PetscObject)ksp), PETSC_ERR_NOT_CONVERGED, "KSPSolve has not converged due to Nan or Inf inner product"); \
498:         { \
499:           PCFailedReason pcreason; \
500:           PetscInt       sendbuf, recvbuf; \
501:           PetscCall(PCGetFailedReasonRank(ksp->pc, &pcreason)); \
502:           sendbuf = (PetscInt)pcreason; \
503:           PetscCall(MPIU_Allreduce(&sendbuf, &recvbuf, 1, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject)ksp))); \
504:           if (recvbuf) { \
505:             PetscCall(PCSetFailedReason(ksp->pc, (PCFailedReason)recvbuf)); \
506:             ksp->reason = KSP_DIVERGED_PC_FAILED; \
507:             PetscCall(VecSetInf(ksp->vec_sol)); \
508:           } else { \
509:             ksp->reason = KSP_DIVERGED_NANORINF; \
510:           } \
511:           PetscFunctionReturn(PETSC_SUCCESS); \
512:         } \
513:       } \
514:     } while (0)

516:   /*MC
517:    KSPCheckNorm - Checks if the result of a norm used by the corresponding `KSP` contains `inf` or `NaN`. These indicate that the previous
518:       application of the preconditioner generated an error. Sets a `KSPConvergedReason` and returns if the `PC` set a `PCFailedReason`.

520:    Collective

522:    Input Parameter:
523: .  ksp - the linear solver `KSP` context.

525:    Output Parameter:
526: .  beta - the result of the norm

528:    Level: developer

530:    Developer Notes:
531:    Used to manage returning from `KSP` solvers whose preconditioners have failed, possibly only a subset of MPI ranks, in some way.

533:    It uses the fact that `KSP` piggy-backs the collectivity of certain error conditions on the results of norms and inner products.

535: .seealso: `PCFailedReason`, `KSPConvergedReason`, `PCGetFailedReasonRank()`, `KSP`, `KSPCreate()`, `KSPSetType()`, `KSP`, `KSPCheckDot()`, `KSPCheckSolve()`
536: M*/
537:   #define KSPCheckNorm(ksp, beta) \
538:     do { \
539:       if (PetscIsInfOrNanReal(beta)) { \
540:         PetscCheck(!ksp->errorifnotconverged, PetscObjectComm((PetscObject)ksp), PETSC_ERR_NOT_CONVERGED, "KSPSolve has not converged due to Nan or Inf norm"); \
541:         { \
542:           PCFailedReason pcreason; \
543:           PetscInt       sendbuf, recvbuf; \
544:           PetscCall(PCGetFailedReasonRank(ksp->pc, &pcreason)); \
545:           sendbuf = (PetscInt)pcreason; \
546:           PetscCall(MPIU_Allreduce(&sendbuf, &recvbuf, 1, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject)ksp))); \
547:           if (recvbuf) { \
548:             PetscCall(PCSetFailedReason(ksp->pc, (PCFailedReason)recvbuf)); \
549:             ksp->reason = KSP_DIVERGED_PC_FAILED; \
550:             PetscCall(VecSetInf(ksp->vec_sol)); \
551:             ksp->rnorm = beta; \
552:           } else { \
553:             PetscCall(PCSetFailedReason(ksp->pc, PC_NOERROR)); \
554:             ksp->reason = KSP_DIVERGED_NANORINF; \
555:             ksp->rnorm  = beta; \
556:           } \
557:           PetscFunctionReturn(PETSC_SUCCESS); \
558:         } \
559:       } \
560:     } while (0)

562: #endif

564: PETSC_INTERN PetscErrorCode KSPMonitorMakeKey_Internal(const char[], PetscViewerType, PetscViewerFormat, char[]);
565: PETSC_INTERN PetscErrorCode KSPMonitorRange_Private(KSP, PetscInt, PetscReal *);