Actual source code: itcreate.c

  1: /*
  2:      The basic KSP routines, Create, View etc. are here.
  3: */
  4: #include <petsc/private/kspimpl.h>

  6: /* Logging support */
  7: PetscClassId  KSP_CLASSID;
  8: PetscClassId  DMKSP_CLASSID;
  9: PetscClassId  KSPGUESS_CLASSID;
 10: PetscLogEvent KSP_GMRESOrthogonalization, KSP_SetUp, KSP_Solve, KSP_SolveTranspose, KSP_MatSolve, KSP_MatSolveTranspose;

 12: /*
 13:    Contains the list of registered KSP routines
 14: */
 15: PetscFunctionList KSPList              = NULL;
 16: PetscBool         KSPRegisterAllCalled = PETSC_FALSE;

 18: /*
 19:    Contains the list of registered KSP monitors
 20: */
 21: PetscFunctionList KSPMonitorList              = NULL;
 22: PetscFunctionList KSPMonitorCreateList        = NULL;
 23: PetscFunctionList KSPMonitorDestroyList       = NULL;
 24: PetscBool         KSPMonitorRegisterAllCalled = PETSC_FALSE;

 26: /*@
 27:   KSPLoad - Loads a `KSP` that has been stored in a `PETSCVIEWERBINARY`  with `KSPView()`.

 29:   Collective

 31:   Input Parameters:
 32: + newdm  - the newly loaded `KSP`, this needs to have been created with `KSPCreate()` or
 33:            some related function before a call to `KSPLoad()`.
 34: - viewer - binary file viewer, obtained from `PetscViewerBinaryOpen()`

 36:   Level: intermediate

 38:   Note:
 39:   The type is determined by the data in the file, any type set into the `KSP` before this call is ignored.

 41: .seealso: [](ch_ksp), `KSP`, `PetscViewerBinaryOpen()`, `KSPView()`, `MatLoad()`, `VecLoad()`
 42: @*/
 43: PetscErrorCode KSPLoad(KSP newdm, PetscViewer viewer)
 44: {
 45:   PetscBool isbinary;
 46:   PetscInt  classid;
 47:   char      type[256];
 48:   PC        pc;

 50:   PetscFunctionBegin;
 53:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
 54:   PetscCheck(isbinary, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid viewer; open viewer with PetscViewerBinaryOpen()");

 56:   PetscCall(PetscViewerBinaryRead(viewer, &classid, 1, NULL, PETSC_INT));
 57:   PetscCheck(classid == KSP_FILE_CLASSID, PetscObjectComm((PetscObject)newdm), PETSC_ERR_ARG_WRONG, "Not KSP next in file");
 58:   PetscCall(PetscViewerBinaryRead(viewer, type, 256, NULL, PETSC_CHAR));
 59:   PetscCall(KSPSetType(newdm, type));
 60:   PetscTryTypeMethod(newdm, load, viewer);
 61:   PetscCall(KSPGetPC(newdm, &pc));
 62:   PetscCall(PCLoad(pc, viewer));
 63:   PetscFunctionReturn(PETSC_SUCCESS);
 64: }

 66: #include <petscdraw.h>
 67: #if defined(PETSC_HAVE_SAWS)
 68: #include <petscviewersaws.h>
 69: #endif
 70: /*@
 71:   KSPView - Prints the various parameters currently set in the `KSP` object. For example, the convergence tolerances and `KSPType`.
 72:   Also views the `PC` and `Mat` contained by the `KSP` with `PCView()` and `MatView()`.

 74:   Collective

 76:   Input Parameters:
 77: + ksp    - the Krylov space context
 78: - viewer - visualization context

 80:   Options Database Key:
 81: . -ksp_view - print the `KSP` data structure at the end of each `KSPSolve()` call

 83:   Level: beginner

 85:   Notes:
 86:   The available visualization contexts include
 87: +     `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
 88: -     `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard
 89:   output where only the first processor opens
 90:   the file.  All other processors send their
 91:   data to the first processor to print.

 93:   The available formats include
 94: +     `PETSC_VIEWER_DEFAULT` - standard output (default)
 95: -     `PETSC_VIEWER_ASCII_INFO_DETAIL` - more verbose output for `PCBJACOBI` and `PCASM`

 97:   The user can open an alternative visualization context with
 98:   `PetscViewerASCIIOpen()` - output to a specified file.

100:   Use `KSPViewFromOptions()` to allow the user to select many different `PetscViewerType` and formats from the options database.

102:   In the debugger you can do call `KSPView(ksp,0)` to display the `KSP`. (The same holds for any PETSc object viewer).

104: .seealso: [](ch_ksp), `KSP`, `PetscViewer`, `PCView()`, `PetscViewerASCIIOpen()`, `KSPViewFromOptions()`
105: @*/
106: PetscErrorCode KSPView(KSP ksp, PetscViewer viewer)
107: {
108:   PetscBool iascii, isbinary, isdraw, isstring;
109: #if defined(PETSC_HAVE_SAWS)
110:   PetscBool issaws;
111: #endif

113:   PetscFunctionBegin;
115:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)ksp), &viewer));
117:   PetscCheckSameComm(ksp, 1, viewer, 2);

119:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
120:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
121:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERDRAW, &isdraw));
122:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSTRING, &isstring));
123: #if defined(PETSC_HAVE_SAWS)
124:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSAWS, &issaws));
125: #endif
126:   if (iascii) {
127:     PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)ksp, viewer));
128:     PetscCall(PetscViewerASCIIPushTab(viewer));
129:     PetscTryTypeMethod(ksp, view, viewer);
130:     PetscCall(PetscViewerASCIIPopTab(viewer));
131:     if (ksp->guess_zero) {
132:       PetscCall(PetscViewerASCIIPrintf(viewer, "  maximum iterations=%" PetscInt_FMT ", initial guess is zero\n", ksp->max_it));
133:     } else {
134:       PetscCall(PetscViewerASCIIPrintf(viewer, "  maximum iterations=%" PetscInt_FMT ", nonzero initial guess\n", ksp->max_it));
135:     }
136:     if (ksp->min_it) PetscCall(PetscViewerASCIIPrintf(viewer, "  minimum iterations=%" PetscInt_FMT "\n", ksp->min_it));
137:     if (ksp->guess_knoll) PetscCall(PetscViewerASCIIPrintf(viewer, "  using preconditioner applied to right-hand side for initial guess\n"));
138:     PetscCall(PetscViewerASCIIPrintf(viewer, "  tolerances: relative=%g, absolute=%g, divergence=%g\n", (double)ksp->rtol, (double)ksp->abstol, (double)ksp->divtol));
139:     if (ksp->pc_side == PC_RIGHT) {
140:       PetscCall(PetscViewerASCIIPrintf(viewer, "  right preconditioning\n"));
141:     } else if (ksp->pc_side == PC_SYMMETRIC) {
142:       PetscCall(PetscViewerASCIIPrintf(viewer, "  symmetric preconditioning\n"));
143:     } else {
144:       PetscCall(PetscViewerASCIIPrintf(viewer, "  left preconditioning\n"));
145:     }
146:     if (ksp->guess) {
147:       PetscCall(PetscViewerASCIIPushTab(viewer));
148:       PetscCall(KSPGuessView(ksp->guess, viewer));
149:       PetscCall(PetscViewerASCIIPopTab(viewer));
150:     }
151:     if (ksp->dscale) PetscCall(PetscViewerASCIIPrintf(viewer, "  diagonally scaled system\n"));
152:     PetscCall(PetscViewerASCIIPrintf(viewer, "  using %s norm type for convergence test\n", KSPNormTypes[ksp->normtype]));
153:   } else if (isbinary) {
154:     PetscInt    classid = KSP_FILE_CLASSID;
155:     MPI_Comm    comm;
156:     PetscMPIInt rank;
157:     char        type[256];

159:     PetscCall(PetscObjectGetComm((PetscObject)ksp, &comm));
160:     PetscCallMPI(MPI_Comm_rank(comm, &rank));
161:     if (rank == 0) {
162:       PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT));
163:       PetscCall(PetscStrncpy(type, ((PetscObject)ksp)->type_name, 256));
164:       PetscCall(PetscViewerBinaryWrite(viewer, type, 256, PETSC_CHAR));
165:     }
166:     PetscTryTypeMethod(ksp, view, viewer);
167:   } else if (isstring) {
168:     const char *type;
169:     PetscCall(KSPGetType(ksp, &type));
170:     PetscCall(PetscViewerStringSPrintf(viewer, " KSPType: %-7.7s", type));
171:     PetscTryTypeMethod(ksp, view, viewer);
172:   } else if (isdraw) {
173:     PetscDraw draw;
174:     char      str[36];
175:     PetscReal x, y, bottom, h;
176:     PetscBool flg;

178:     PetscCall(PetscViewerDrawGetDraw(viewer, 0, &draw));
179:     PetscCall(PetscDrawGetCurrentPoint(draw, &x, &y));
180:     PetscCall(PetscObjectTypeCompare((PetscObject)ksp, KSPPREONLY, &flg));
181:     if (!flg) {
182:       PetscCall(PetscStrncpy(str, "KSP: ", sizeof(str)));
183:       PetscCall(PetscStrlcat(str, ((PetscObject)ksp)->type_name, sizeof(str)));
184:       PetscCall(PetscDrawStringBoxed(draw, x, y, PETSC_DRAW_RED, PETSC_DRAW_BLACK, str, NULL, &h));
185:       bottom = y - h;
186:     } else {
187:       bottom = y;
188:     }
189:     PetscCall(PetscDrawPushCurrentPoint(draw, x, bottom));
190: #if defined(PETSC_HAVE_SAWS)
191:   } else if (issaws) {
192:     PetscMPIInt rank;
193:     const char *name;

195:     PetscCall(PetscObjectGetName((PetscObject)ksp, &name));
196:     PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
197:     if (!((PetscObject)ksp)->amsmem && rank == 0) {
198:       char dir[1024];

200:       PetscCall(PetscObjectViewSAWs((PetscObject)ksp, viewer));
201:       PetscCall(PetscSNPrintf(dir, 1024, "/PETSc/Objects/%s/its", name));
202:       PetscCallSAWs(SAWs_Register, (dir, &ksp->its, 1, SAWs_READ, SAWs_INT));
203:       if (!ksp->res_hist) PetscCall(KSPSetResidualHistory(ksp, NULL, PETSC_DECIDE, PETSC_TRUE));
204:       PetscCall(PetscSNPrintf(dir, 1024, "/PETSc/Objects/%s/res_hist", name));
205:       PetscCallSAWs(SAWs_Register, (dir, ksp->res_hist, 10, SAWs_READ, SAWs_DOUBLE));
206:     }
207: #endif
208:   } else PetscTryTypeMethod(ksp, view, viewer);
209:   if (ksp->pc) PetscCall(PCView(ksp->pc, viewer));
210:   if (isdraw) {
211:     PetscDraw draw;
212:     PetscCall(PetscViewerDrawGetDraw(viewer, 0, &draw));
213:     PetscCall(PetscDrawPopCurrentPoint(draw));
214:   }
215:   PetscFunctionReturn(PETSC_SUCCESS);
216: }

218: /*@
219:   KSPViewFromOptions - View (print) a `KSP` object based on values in the options database. Also views the `PC` and `Mat` contained by the `KSP`
220:   with `PCView()` and `MatView()`.

222:   Collective

224:   Input Parameters:
225: + A    - Krylov solver context
226: . obj  - Optional object that provides the options prefix used to query the options database
227: - name - command line option

229:   Level: intermediate

231: .seealso: [](ch_ksp), `KSP`, `KSPView()`, `PetscObjectViewFromOptions()`, `KSPCreate()`
232: @*/
233: PetscErrorCode KSPViewFromOptions(KSP A, PetscObject obj, const char name[])
234: {
235:   PetscFunctionBegin;
237:   PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
238:   PetscFunctionReturn(PETSC_SUCCESS);
239: }

241: /*@
242:   KSPSetNormType - Sets the type of residual norm that is used for convergence testing in `KSPSolve()` for the given `KSP` context

244:   Logically Collective

246:   Input Parameters:
247: + ksp      - Krylov solver context
248: - normtype - one of
249: .vb
250:    KSP_NORM_NONE             - skips computing the norm, this should generally only be used if you are using
251:                                the Krylov method as a smoother with a fixed small number of iterations.
252:                                Implicitly sets `KSPConvergedSkip()` as the `KSP` convergence test.
253:                                Note that certain algorithms such as `KSPGMRES` ALWAYS require the norm calculation,
254:                                for these methods the norms are still computed, they are just not used in
255:                                the convergence test.
256:    KSP_NORM_PRECONDITIONED   - the default for left-preconditioned solves, uses the 2-norm
257:                                of the preconditioned residual  $B^{-1}(b - A x)$.
258:    KSP_NORM_UNPRECONDITIONED - uses the 2-norm of the true $b - Ax$ residual.
259:    KSP_NORM_NATURAL          - uses the $A$ norm of the true $b - Ax$ residual; supported by `KSPCG`, `KSPCR`, `KSPCGNE`, `KSPCGS`
260: .ve

262:   Options Database Key:
263: . -ksp_norm_type <none,preconditioned,unpreconditioned,natural> - set `KSP` norm type

265:   Level: advanced

267:   Notes:
268:   The norm is always of the equations residual $\| b - A x^n \|$  (or an approximation to that norm), they are never a norm of the error in the equation.

270:   Not all combinations of preconditioner side (see `KSPSetPCSide()`) and norm types are supported by all Krylov methods.
271:   If only one is set, PETSc tries to automatically change the other to find a compatible pair.  If no such combination
272:   is supported, PETSc will generate an error.

274:   Developer Note:
275:   Supported combinations of norm and preconditioner side are set using `KSPSetSupportedNorm()` for each `KSPType`.

277: .seealso: [](ch_ksp), `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSPConvergedSkip()`, `KSPSetCheckNormIteration()`, `KSPSetPCSide()`, `KSPGetPCSide()`, `KSPNormType`
278: @*/
279: PetscErrorCode KSPSetNormType(KSP ksp, KSPNormType normtype)
280: {
281:   PetscFunctionBegin;
284:   ksp->normtype = ksp->normtype_set = normtype;
285:   PetscFunctionReturn(PETSC_SUCCESS);
286: }

288: /*@
289:   KSPSetCheckNormIteration - Sets the first iteration at which the norm of the residual will be
290:   computed and used in the convergence test of `KSPSolve()` for the given `KSP` context

292:   Logically Collective

294:   Input Parameters:
295: + ksp - Krylov solver context
296: - it  - use -1 to check at all iterations

298:   Level: advanced

300:   Notes:
301:   Currently only works with `KSPCG`, `KSPBCGS` and `KSPIBCGS`

303:   Use `KSPSetNormType`(ksp,`KSP_NORM_NONE`) to never check the norm

305:   On steps where the norm is not computed, the previous norm is still in the variable, so if you run with, for example,
306:   `-ksp_monitor` the residual norm will appear to be unchanged for several iterations (though it is not really unchanged).

308:   Certain methods such as `KSPGMRES` always compute the residual norm, this routine will not change that computation, but it will
309:   prevent the computed norm from being checked.

311: .seealso: [](ch_ksp), `KSP`, `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSPConvergedSkip()`, `KSPSetNormType()`, `KSPSetLagNorm()`
312: @*/
313: PetscErrorCode KSPSetCheckNormIteration(KSP ksp, PetscInt it)
314: {
315:   PetscFunctionBegin;
318:   ksp->chknorm = it;
319:   PetscFunctionReturn(PETSC_SUCCESS);
320: }

322: /*@
323:   KSPSetLagNorm - Lags the residual norm calculation so that it is computed as part of the `MPI_Allreduce()` used for
324:   computing the inner products needed for the next iteration.

326:   Logically Collective

328:   Input Parameters:
329: + ksp - Krylov solver context
330: - flg - `PETSC_TRUE` or `PETSC_FALSE`

332:   Options Database Key:
333: . -ksp_lag_norm - lag the calculated residual norm

335:   Level: advanced

337:   Notes:
338:   Currently only works with `KSPIBCGS`.

340:   This can reduce communication costs at the expense of doing
341:   one additional iteration because the norm used in the convergence test of `KSPSolve()` is one iteration behind the actual
342:   current residual norm (which has not yet been computed due to the lag).

344:   Use `KSPSetNormType`(ksp,`KSP_NORM_NONE`) to never check the norm

346:   If you lag the norm and run with, for example, `-ksp_monitor`, the residual norm reported will be the lagged one.

348:   `KSPSetCheckNormIteration()` is an alternative way of avoiding the expense of computing the residual norm at each iteration.

350: .seealso: [](ch_ksp), `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSPConvergedSkip()`, `KSPSetNormType()`, `KSPSetCheckNormIteration()`
351: @*/
352: PetscErrorCode KSPSetLagNorm(KSP ksp, PetscBool flg)
353: {
354:   PetscFunctionBegin;
357:   ksp->lagnorm = flg;
358:   PetscFunctionReturn(PETSC_SUCCESS);
359: }

361: /*@
362:   KSPSetSupportedNorm - Sets a norm and preconditioner side supported by a `KSPType`

364:   Logically Collective

366:   Input Parameters:
367: + ksp      - Krylov method
368: . normtype - supported norm type of the type `KSPNormType`
369: . pcside   - preconditioner side, of the type `PCSide` that can be used with this `KSPNormType`
370: - priority - positive integer preference for this combination; larger values have higher priority

372:   Level: developer

374:   Notes:
375:   This function should be called from the implementation files `KSPCreate_XXX()` to declare
376:   which norms and preconditioner sides are supported. Users should not call this
377:   function.

379:   This function can be called multiple times for each combination of `KSPNormType` and `PCSide`
380:   the `KSPType` supports

382: .seealso: [](ch_ksp), `KSP`, `KSPNormType`, `PCSide`, `KSPSetNormType()`, `KSPSetPCSide()`
383: @*/
384: PetscErrorCode KSPSetSupportedNorm(KSP ksp, KSPNormType normtype, PCSide pcside, PetscInt priority)
385: {
386:   PetscFunctionBegin;
388:   ksp->normsupporttable[normtype][pcside] = priority;
389:   PetscFunctionReturn(PETSC_SUCCESS);
390: }

392: static PetscErrorCode KSPNormSupportTableReset_Private(KSP ksp)
393: {
394:   PetscFunctionBegin;
395:   PetscCall(PetscMemzero(ksp->normsupporttable, sizeof(ksp->normsupporttable)));
396:   ksp->pc_side  = ksp->pc_side_set;
397:   ksp->normtype = ksp->normtype_set;
398:   PetscFunctionReturn(PETSC_SUCCESS);
399: }

401: PetscErrorCode KSPSetUpNorms_Private(KSP ksp, PetscBool errorifnotsupported, KSPNormType *normtype, PCSide *pcside)
402: {
403:   PetscInt i, j, best, ibest = 0, jbest = 0;

405:   PetscFunctionBegin;
406:   best = 0;
407:   for (i = 0; i < KSP_NORM_MAX; i++) {
408:     for (j = 0; j < PC_SIDE_MAX; j++) {
409:       if ((ksp->normtype == KSP_NORM_DEFAULT || ksp->normtype == i) && (ksp->pc_side == PC_SIDE_DEFAULT || ksp->pc_side == j) && (ksp->normsupporttable[i][j] > best)) {
410:         best  = ksp->normsupporttable[i][j];
411:         ibest = i;
412:         jbest = j;
413:       }
414:     }
415:   }
416:   if (best < 1 && errorifnotsupported) {
417:     PetscCheck(ksp->normtype != KSP_NORM_DEFAULT || ksp->pc_side != PC_SIDE_DEFAULT, PetscObjectComm((PetscObject)ksp), PETSC_ERR_PLIB, "The %s KSP implementation did not call KSPSetSupportedNorm()", ((PetscObject)ksp)->type_name);
418:     PetscCheck(ksp->normtype != KSP_NORM_DEFAULT, PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "KSP %s does not support preconditioner side %s", ((PetscObject)ksp)->type_name, PCSides[ksp->pc_side]);
419:     PetscCheck(ksp->pc_side != PC_SIDE_DEFAULT, PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "KSP %s does not support norm type %s", ((PetscObject)ksp)->type_name, KSPNormTypes[ksp->normtype]);
420:     SETERRQ(PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "KSP %s does not support norm type %s with preconditioner side %s", ((PetscObject)ksp)->type_name, KSPNormTypes[ksp->normtype], PCSides[ksp->pc_side]);
421:   }
422:   if (normtype) *normtype = (KSPNormType)ibest;
423:   if (pcside) *pcside = (PCSide)jbest;
424:   PetscFunctionReturn(PETSC_SUCCESS);
425: }

427: /*@
428:   KSPGetNormType - Gets the `KSPNormType` that is used for convergence testing during `KSPSolve()` for this `KSP` context

430:   Not Collective

432:   Input Parameter:
433: . ksp - Krylov solver context

435:   Output Parameter:
436: . normtype - the `KSPNormType` that is used for convergence testing

438:   Level: advanced

440: .seealso: [](ch_ksp), `KSPNormType`, `KSPSetNormType()`, `KSPConvergedSkip()`
441: @*/
442: PetscErrorCode KSPGetNormType(KSP ksp, KSPNormType *normtype)
443: {
444:   PetscFunctionBegin;
446:   PetscAssertPointer(normtype, 2);
447:   PetscCall(KSPSetUpNorms_Private(ksp, PETSC_TRUE, &ksp->normtype, &ksp->pc_side));
448:   *normtype = ksp->normtype;
449:   PetscFunctionReturn(PETSC_SUCCESS);
450: }

452: #if defined(PETSC_HAVE_SAWS)
453: #include <petscviewersaws.h>
454: #endif

456: /*@
457:   KSPSetOperators - Sets the matrix associated with the linear system
458:   and a (possibly) different one from which the preconditioner will be built into the `KSP` context. The matrix will then be used during `KSPSolve()`

460:   Collective

462:   Input Parameters:
463: + ksp  - the `KSP` context
464: . Amat - the matrix that defines the linear system
465: - Pmat - the matrix to be used in constructing the preconditioner, usually the same as `Amat`.

467:   Level: beginner

469:   Notes:
470:   If you know the operator `Amat` has a null space you can use `MatSetNullSpace()` and `MatSetTransposeNullSpace()` to supply the null
471:   space to `Amat` and the `KSP` solvers will automatically use that null space as needed during the solution process.

473:   All future calls to `KSPSetOperators()` must use the same size matrices, unless `KSPReset()` is called!

475:   Passing a `NULL` for `Amat` or `Pmat` removes the matrix that is currently being used from the `KSP` context.

477:   If you wish to replace either `Amat` or `Pmat` but leave the other one untouched then
478:   first call `KSPGetOperators()` to get the one you wish to keep, call `PetscObjectReference()`
479:   on it and then pass it back in your call to `KSPSetOperators()`.

481:   Developer Notes:
482:   If the operators have NOT been set with `KSPSetOperators()` then the operators
483:   are created in the `PC` and returned to the user. In this case, if both operators
484:   mat and pmat are requested, two DIFFERENT operators will be returned. If
485:   only one is requested both operators in the `PC` will be the same (i.e. as
486:   if one had called `KSPSetOperators()` with the same argument for both `Mat`s).
487:   The user must set the sizes of the returned matrices and their type etc just
488:   as if the user created them with `MatCreate()`. For example,

490: .vb
491:          KSPGetOperators(ksp/pc,&mat,NULL); is equivalent to
492:            set size, type, etc of mat

494:          MatCreate(comm,&mat);
495:          KSP/PCSetOperators(ksp/pc,mat,mat);
496:          PetscObjectDereference((PetscObject)mat);
497:            set size, type, etc of mat

499:      and

501:          KSP/PCGetOperators(ksp/pc,&mat,&pmat); is equivalent to
502:            set size, type, etc of mat and pmat

504:          MatCreate(comm,&mat);
505:          MatCreate(comm,&pmat);
506:          KSP/PCSetOperators(ksp/pc,mat,pmat);
507:          PetscObjectDereference((PetscObject)mat);
508:          PetscObjectDereference((PetscObject)pmat);
509:            set size, type, etc of mat and pmat
510: .ve

512:   The rationale for this support is so that when creating a `TS`, `SNES`, or `KSP` the hierarchy
513:   of underlying objects (i.e. `SNES`, `KSP`, `PC`, `Mat`) and their lifespans can be completely
514:   managed by the top most level object (i.e. the `TS`, `SNES`, or `KSP`). Another way to look
515:   at this is when you create a `SNES` you do not NEED to create a `KSP` and attach it to
516:   the `SNES` object (the `SNES` object manages it for you). Similarly when you create a `KSP`
517:   you do not need to attach a `PC` to it (the `KSP` object manages the `PC` object for you).
518:   Thus, why should YOU have to create the `Mat` and attach it to the `SNES`/`KSP`/`PC`, when
519:   it can be created for you?

521: .seealso: [](ch_ksp), `KSP`, `Mat`, `KSPSolve()`, `KSPGetPC()`, `PCGetOperators()`, `PCSetOperators()`, `KSPGetOperators()`, `KSPSetComputeOperators()`, `KSPSetComputeInitialGuess()`, `KSPSetComputeRHS()`
522: @*/
523: PetscErrorCode KSPSetOperators(KSP ksp, Mat Amat, Mat Pmat)
524: {
525:   PetscFunctionBegin;
529:   if (Amat) PetscCheckSameComm(ksp, 1, Amat, 2);
530:   if (Pmat) PetscCheckSameComm(ksp, 1, Pmat, 3);
531:   if (!ksp->pc) PetscCall(KSPGetPC(ksp, &ksp->pc));
532:   PetscCall(PCSetOperators(ksp->pc, Amat, Pmat));
533:   if (ksp->setupstage == KSP_SETUP_NEWRHS) ksp->setupstage = KSP_SETUP_NEWMATRIX; /* so that next solve call will call PCSetUp() on new matrix */
534:   PetscFunctionReturn(PETSC_SUCCESS);
535: }

537: /*@
538:   KSPGetOperators - Gets the matrix associated with the linear system
539:   and a (possibly) different one used to construct the preconditioner from the `KSP` context

541:   Collective

543:   Input Parameter:
544: . ksp - the `KSP` context

546:   Output Parameters:
547: + Amat - the matrix that defines the linear system
548: - Pmat - the matrix to be used in constructing the preconditioner, usually the same as `Amat`.

550:   Level: intermediate

552:   Notes:
553:   If `KSPSetOperators()` has not been called then the `KSP` object will attempt to automatically create the matrix `Amat` and return it

555:   Use `KSPGetOperatorsSet()` to determine if matrices have been provided.

557:   DOES NOT increase the reference counts of the matrix, so you should NOT destroy them.

559: .seealso: [](ch_ksp), `KSP`, `KSPSolve()`, `KSPGetPC()`, `PCSetOperators()`, `KSPSetOperators()`, `KSPGetOperatorsSet()`
560: @*/
561: PetscErrorCode KSPGetOperators(KSP ksp, Mat *Amat, Mat *Pmat)
562: {
563:   PetscFunctionBegin;
565:   if (!ksp->pc) PetscCall(KSPGetPC(ksp, &ksp->pc));
566:   PetscCall(PCGetOperators(ksp->pc, Amat, Pmat));
567:   PetscFunctionReturn(PETSC_SUCCESS);
568: }

570: /*@
571:   KSPGetOperatorsSet - Determines if the matrix associated with the linear system and
572:   possibly a different one from which the preconditioner will be built have been set in the `KSP` with `KSPSetOperators()`

574:   Not Collective, though the results on all processes will be the same

576:   Input Parameter:
577: . ksp - the `KSP` context

579:   Output Parameters:
580: + mat  - the matrix associated with the linear system was set
581: - pmat - matrix from which the preconditioner will be built, usually the same as `mat` was set

583:   Level: intermediate

585:   Note:
586:   This routine exists because if you call `KSPGetOperators()` on a `KSP` that does not yet have operators they are
587:   automatically created in the call.

589: .seealso: [](ch_ksp), `KSP`, `PCSetOperators()`, `KSPGetOperators()`, `KSPSetOperators()`, `PCGetOperators()`, `PCGetOperatorsSet()`
590: @*/
591: PetscErrorCode KSPGetOperatorsSet(KSP ksp, PetscBool *mat, PetscBool *pmat)
592: {
593:   PetscFunctionBegin;
595:   if (!ksp->pc) PetscCall(KSPGetPC(ksp, &ksp->pc));
596:   PetscCall(PCGetOperatorsSet(ksp->pc, mat, pmat));
597:   PetscFunctionReturn(PETSC_SUCCESS);
598: }

600: /*@C
601:   KSPSetPreSolve - Sets a function that is called at the beginning of each `KSPSolve()`. Used in conjunction with `KSPSetPostSolve()`.

603:   Logically Collective

605:   Input Parameters:
606: + ksp      - the solver object
607: . presolve - the function to call before the solve
608: - ctx      - an optional context needed by the function

610:   Calling sequence of `presolve`:
611: + ksp - the `KSP` context
612: . rhs - the right-hand side vector
613: . x   - the solution vector
614: - ctx - optional user-provided context

616:   Level: developer

618:   Notes:
619:   The function provided here `presolve` is used to modify the right hand side, and possibly the matrix, of the linear system to be solved.
620:   The function provided with `KSPSetPostSolve()` then modifies the resulting solution of that linear system to obtain the correct solution
621:   to the initial linear system.

623:   The functions `PCPreSolve()` and `PCPostSolve()` provide a similar functionality and are used, for example with `PCEISENSTAT`.

625: .seealso: [](ch_ksp), `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPSetPostSolve()`, `PCEISENSTAT`, `PCPreSolve()`, `PCPostSolve()`
626: @*/
627: PetscErrorCode KSPSetPreSolve(KSP ksp, PetscErrorCode (*presolve)(KSP ksp, Vec rhs, Vec x, void *ctx), void *ctx)
628: {
629:   PetscFunctionBegin;
631:   ksp->presolve = presolve;
632:   ksp->prectx   = ctx;
633:   PetscFunctionReturn(PETSC_SUCCESS);
634: }

636: /*@C
637:   KSPSetPostSolve - Sets a function that is called at the end of each `KSPSolve()` (whether it converges or not). Used in conjunction with `KSPSetPreSolve()`.

639:   Logically Collective

641:   Input Parameters:
642: + ksp       - the solver object
643: . postsolve - the function to call after the solve
644: - ctx       - an optional context needed by the function

646:   Calling sequence of `postsolve`:
647: + ksp - the `KSP` context
648: . rhs - the right-hand side vector
649: . x   - the solution vector
650: - ctx - optional user-provided context

652:   Level: developer

654: .seealso: [](ch_ksp), `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPSetPreSolve()`, `PCEISENSTAT`
655: @*/
656: PetscErrorCode KSPSetPostSolve(KSP ksp, PetscErrorCode (*postsolve)(KSP ksp, Vec rhs, Vec x, void *ctx), void *ctx)
657: {
658:   PetscFunctionBegin;
660:   ksp->postsolve = postsolve;
661:   ksp->postctx   = ctx;
662:   PetscFunctionReturn(PETSC_SUCCESS);
663: }

665: /*@
666:   KSPSetNestLevel - sets the amount of nesting the `KSP` has. That is the number of levels of `KSP` above this `KSP` in a linear solve.

668:   Collective

670:   Input Parameters:
671: + ksp   - the `KSP`
672: - level - the nest level

674:   Level: developer

676:   Note:
677:   For example, the `KSP` in each block of a `KSPBJACOBI` has a level of 1, while the outer `KSP` has a level of 0.

679: .seealso: [](ch_ksp), `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPGMRES`, `KSPType`, `KSPGetNestLevel()`, `PCSetKSPNestLevel()`, `PCGetKSPNestLevel()`
680: @*/
681: PetscErrorCode KSPSetNestLevel(KSP ksp, PetscInt level)
682: {
683:   PetscFunctionBegin;
686:   ksp->nestlevel = level;
687:   PetscFunctionReturn(PETSC_SUCCESS);
688: }

690: /*@
691:   KSPGetNestLevel - gets the amount of nesting the `KSP` has

693:   Not Collective

695:   Input Parameter:
696: . ksp - the `KSP`

698:   Output Parameter:
699: . level - the nest level

701:   Level: developer

703: .seealso: [](ch_ksp), `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPGMRES`, `KSPType`, `KSPSetNestLevel()`, `PCSetKSPNestLevel()`, `PCGetKSPNestLevel()`
704: @*/
705: PetscErrorCode KSPGetNestLevel(KSP ksp, PetscInt *level)
706: {
707:   PetscFunctionBegin;
709:   PetscAssertPointer(level, 2);
710:   *level = ksp->nestlevel;
711:   PetscFunctionReturn(PETSC_SUCCESS);
712: }

714: /*@
715:   KSPCreate - Creates the `KSP` context. This `KSP` context is used in PETSc to solve linear systems with `KSPSolve()`

717:   Collective

719:   Input Parameter:
720: . comm - MPI communicator

722:   Output Parameter:
723: . inksp - location to put the `KSP` context

725:   Level: beginner

727:   Note:
728:   The default `KSPType` is `KSPGMRES` with a restart of 30, using modified Gram-Schmidt orthogonalization. The `KSPType` may be
729:   changed with `KSPSetType()`

731: .seealso: [](ch_ksp), `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPGMRES`, `KSPType`, `KSPSetType()`
732: @*/
733: PetscErrorCode KSPCreate(MPI_Comm comm, KSP *inksp)
734: {
735:   KSP   ksp;
736:   void *ctx;

738:   PetscFunctionBegin;
739:   PetscAssertPointer(inksp, 2);
740:   PetscCall(KSPInitializePackage());

742:   PetscCall(PetscHeaderCreate(ksp, KSP_CLASSID, "KSP", "Krylov Method", "KSP", comm, KSPDestroy, KSPView));
743:   ksp->default_max_it = ksp->max_it = 10000;
744:   ksp->pc_side = ksp->pc_side_set = PC_SIDE_DEFAULT;

746:   ksp->default_rtol = ksp->rtol = 1.e-5;
747:   ksp->default_abstol = ksp->abstol = PetscDefined(USE_REAL_SINGLE) ? 1.e-25 : 1.e-50;
748:   ksp->default_divtol = ksp->divtol = 1.e4;

750:   ksp->chknorm  = -1;
751:   ksp->normtype = ksp->normtype_set = KSP_NORM_DEFAULT;
752:   ksp->rnorm                        = 0.0;
753:   ksp->its                          = 0;
754:   ksp->guess_zero                   = PETSC_TRUE;
755:   ksp->calc_sings                   = PETSC_FALSE;
756:   ksp->res_hist                     = NULL;
757:   ksp->res_hist_alloc               = NULL;
758:   ksp->res_hist_len                 = 0;
759:   ksp->res_hist_max                 = 0;
760:   ksp->res_hist_reset               = PETSC_TRUE;
761:   ksp->err_hist                     = NULL;
762:   ksp->err_hist_alloc               = NULL;
763:   ksp->err_hist_len                 = 0;
764:   ksp->err_hist_max                 = 0;
765:   ksp->err_hist_reset               = PETSC_TRUE;
766:   ksp->numbermonitors               = 0;
767:   ksp->numberreasonviews            = 0;
768:   ksp->setfromoptionscalled         = 0;
769:   ksp->nmax                         = PETSC_DECIDE;

771:   PetscCall(KSPConvergedDefaultCreate(&ctx));
772:   PetscCall(KSPSetConvergenceTest(ksp, KSPConvergedDefault, ctx, KSPConvergedDefaultDestroy));
773:   ksp->ops->buildsolution = KSPBuildSolutionDefault;
774:   ksp->ops->buildresidual = KSPBuildResidualDefault;

776:   ksp->vec_sol    = NULL;
777:   ksp->vec_rhs    = NULL;
778:   ksp->pc         = NULL;
779:   ksp->data       = NULL;
780:   ksp->nwork      = 0;
781:   ksp->work       = NULL;
782:   ksp->reason     = KSP_CONVERGED_ITERATING;
783:   ksp->setupstage = KSP_SETUP_NEW;

785:   PetscCall(KSPNormSupportTableReset_Private(ksp));

787:   *inksp = ksp;
788:   PetscFunctionReturn(PETSC_SUCCESS);
789: }

791: /*@
792:   KSPSetType - Sets the algorithm/method to be used to solve the linear system with the given `KSP`

794:   Logically Collective

796:   Input Parameters:
797: + ksp  - the Krylov space context
798: - type - a known method

800:   Options Database Key:
801: . -ksp_type  <method> - Sets the method; see `KSPGType` or use `-help` for a list  of available methods (for instance, cg or gmres)

803:   Level: intermediate

805:   Notes:
806:   See `KSPType` for available methods (for instance, `KSPCG` or `KSPGMRES`).

808:   Normally, it is best to use the `KSPSetFromOptions()` command and
809:   then set the `KSP` type from the options database rather than by using
810:   this routine.  Using the options database provides the user with
811:   maximum flexibility in evaluating the many different Krylov methods.
812:   The `KSPSetType()` routine is provided for those situations where it
813:   is necessary to set the iterative solver independently of the command
814:   line or options database.  This might be the case, for example, when
815:   the choice of iterative solver changes during the execution of the
816:   program, and the user's application is taking responsibility for
817:   choosing the appropriate method.  In other words, this routine is
818:   not for beginners.

820:   Developer Note:
821:   `KSPRegister()` is used to add Krylov types to `KSPList` from which they are accessed by `KSPSetType()`.

823: .seealso: [](ch_ksp), `PCSetType()`, `KSPType`, `KSPRegister()`, `KSPCreate()`, `KSP`
824: @*/
825: PetscErrorCode KSPSetType(KSP ksp, KSPType type)
826: {
827:   PetscBool match;
828:   PetscErrorCode (*r)(KSP);

830:   PetscFunctionBegin;
832:   PetscAssertPointer(type, 2);

834:   PetscCall(PetscObjectTypeCompare((PetscObject)ksp, type, &match));
835:   if (match) PetscFunctionReturn(PETSC_SUCCESS);

837:   PetscCall(PetscFunctionListFind(KSPList, type, &r));
838:   PetscCheck(r, PetscObjectComm((PetscObject)ksp), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unable to find requested KSP type %s", type);
839:   /* Destroy the previous private KSP context */
840:   PetscTryTypeMethod(ksp, destroy);

842:   /* Reinitialize function pointers in KSPOps structure */
843:   PetscCall(PetscMemzero(ksp->ops, sizeof(struct _KSPOps)));
844:   ksp->ops->buildsolution = KSPBuildSolutionDefault;
845:   ksp->ops->buildresidual = KSPBuildResidualDefault;
846:   PetscCall(KSPNormSupportTableReset_Private(ksp));
847:   ksp->converged_neg_curve = PETSC_FALSE; // restore default
848:   ksp->setupnewmatrix      = PETSC_FALSE; // restore default (setup not called in case of new matrix)
849:   /* Call the KSPCreate_XXX routine for this particular Krylov solver */
850:   ksp->setupstage     = KSP_SETUP_NEW;
851:   ksp->guess_not_read = PETSC_FALSE; // restore default
852:   PetscCall((*r)(ksp));
853:   PetscCall(PetscObjectChangeTypeName((PetscObject)ksp, type));
854:   PetscFunctionReturn(PETSC_SUCCESS);
855: }

857: /*@
858:   KSPGetType - Gets the `KSP` type as a string from the `KSP` object.

860:   Not Collective

862:   Input Parameter:
863: . ksp - Krylov context

865:   Output Parameter:
866: . type - name of the `KSP` method

868:   Level: intermediate

870: .seealso: [](ch_ksp), `KSPType`, `KSP`, `KSPSetType()`
871: @*/
872: PetscErrorCode KSPGetType(KSP ksp, KSPType *type)
873: {
874:   PetscFunctionBegin;
876:   PetscAssertPointer(type, 2);
877:   *type = ((PetscObject)ksp)->type_name;
878:   PetscFunctionReturn(PETSC_SUCCESS);
879: }

881: /*@C
882:   KSPRegister -  Adds a method, `KSPType`, to the Krylov subspace solver package.

884:   Not Collective, No Fortran Support

886:   Input Parameters:
887: + sname    - name of a new user-defined solver
888: - function - routine to create method

890:   Level: advanced

892:   Note:
893:   `KSPRegister()` may be called multiple times to add several user-defined solvers.

895:   Example Usage:
896: .vb
897:    KSPRegister("my_solver", MySolverCreate);
898: .ve

900:   Then, your solver can be chosen with the procedural interface via
901: .vb
902:   KSPSetType(ksp, "my_solver")
903: .ve
904:   or at runtime via the option `-ksp_type my_solver`

906: .seealso: [](ch_ksp), `KSP`, `KSPType`, `KSPSetType`, `KSPRegisterAll()`
907: @*/
908: PetscErrorCode KSPRegister(const char sname[], PetscErrorCode (*function)(KSP))
909: {
910:   PetscFunctionBegin;
911:   PetscCall(KSPInitializePackage());
912:   PetscCall(PetscFunctionListAdd(&KSPList, sname, function));
913:   PetscFunctionReturn(PETSC_SUCCESS);
914: }

916: PetscErrorCode KSPMonitorMakeKey_Internal(const char name[], PetscViewerType vtype, PetscViewerFormat format, char key[])
917: {
918:   PetscFunctionBegin;
919:   PetscCall(PetscStrncpy(key, name, PETSC_MAX_PATH_LEN));
920:   PetscCall(PetscStrlcat(key, ":", PETSC_MAX_PATH_LEN));
921:   PetscCall(PetscStrlcat(key, vtype, PETSC_MAX_PATH_LEN));
922:   PetscCall(PetscStrlcat(key, ":", PETSC_MAX_PATH_LEN));
923:   PetscCall(PetscStrlcat(key, PetscViewerFormats[format], PETSC_MAX_PATH_LEN));
924:   PetscFunctionReturn(PETSC_SUCCESS);
925: }

927: /*@C
928:   KSPMonitorRegister -  Registers a Krylov subspace solver monitor routine that may be accessed with `KSPMonitorSetFromOptions()`

930:   Not Collective

932:   Input Parameters:
933: + name    - name of a new monitor routine
934: . vtype   - A `PetscViewerType` for the output
935: . format  - A `PetscViewerFormat` for the output
936: . monitor - Monitor routine
937: . create  - Creation routine, or `NULL`
938: - destroy - Destruction routine, or `NULL`

940:   Level: advanced

942:   Notes:
943:   `KSPMonitorRegister()` may be called multiple times to add several user-defined monitors.

945:   The calling sequence for the given function matches the calling sequence used by functions passed to `KSPMonitorSet()` with the additional
946:   requirement that its final argument be a `PetscViewerAndFormat`.

948:   Example Usage:
949: .vb
950:   KSPMonitorRegister("my_monitor", PETSCVIEWERASCII, PETSC_VIEWER_ASCII_INFO_DETAIL, MyMonitor, NULL, NULL);
951: .ve

953:   Then, your monitor can be chosen with the procedural interface via
954: .vb
955:   KSPMonitorSetFromOptions(ksp, "-ksp_monitor_my_monitor", "my_monitor", NULL)
956: .ve
957:   or at runtime via the option `-ksp_monitor_my_monitor`

959: .seealso: [](ch_ksp), `KSP`, `KSPMonitorSet()`, `KSPMonitorRegisterAll()`, `KSPMonitorSetFromOptions()`
960: @*/
961: PetscErrorCode KSPMonitorRegister(const char name[], PetscViewerType vtype, PetscViewerFormat format, PetscErrorCode (*monitor)(KSP, PetscInt, PetscReal, PetscViewerAndFormat *), PetscErrorCode (*create)(PetscViewer, PetscViewerFormat, void *, PetscViewerAndFormat **), PetscErrorCode (*destroy)(PetscViewerAndFormat **))
962: {
963:   char key[PETSC_MAX_PATH_LEN];

965:   PetscFunctionBegin;
966:   PetscCall(KSPInitializePackage());
967:   PetscCall(KSPMonitorMakeKey_Internal(name, vtype, format, key));
968:   PetscCall(PetscFunctionListAdd(&KSPMonitorList, key, monitor));
969:   if (create) PetscCall(PetscFunctionListAdd(&KSPMonitorCreateList, key, create));
970:   if (destroy) PetscCall(PetscFunctionListAdd(&KSPMonitorDestroyList, key, destroy));
971:   PetscFunctionReturn(PETSC_SUCCESS);
972: }