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: /*@C
 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: `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: /*@C
 71:    KSPView - Prints the `KSP` data structure.

 73:    Collective

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

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

 82:    Level: beginner

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

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

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

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

101: .seealso: `KSP`, `PetscViewer`, `PCView()`, `PetscViewerASCIIOpen()`
102: @*/
103: PetscErrorCode KSPView(KSP ksp, PetscViewer viewer)
104: {
105:   PetscBool iascii, isbinary, isdraw, isstring;
106: #if defined(PETSC_HAVE_SAWS)
107:   PetscBool issaws;
108: #endif

110:   PetscFunctionBegin;
112:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)ksp), &viewer));
114:   PetscCheckSameComm(ksp, 1, viewer, 2);

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

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

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

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

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

215: /*@C
216:    KSPViewFromOptions - View a `KSP` object based on values in the options database

218:    Collective

220:    Input Parameters:
221: +  A - Krylov solver context
222: .  obj - Optional object
223: -  name - command line option

225:    Level: intermediate

227: .seealso: `KSP`, `KSPView`, `PetscObjectViewFromOptions()`, `KSPCreate()`
228: @*/
229: PetscErrorCode KSPViewFromOptions(KSP A, PetscObject obj, const char name[])
230: {
231:   PetscFunctionBegin;
233:   PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
234:   PetscFunctionReturn(PETSC_SUCCESS);
235: }

237: /*@
238:    KSPSetNormType - Sets the norm that is used for convergence testing.

240:    Logically Collective

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

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

261:    Level: advanced

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

268:    Developer Note:
269:    Supported combinations of norm and preconditioner side are set using `KSPSetSupportedNorm()`.

271: .seealso: `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSPConvergedSkip()`, `KSPSetCheckNormIteration()`, `KSPSetPCSide()`, `KSPGetPCSide()`, `KSPNormType`
272: @*/
273: PetscErrorCode KSPSetNormType(KSP ksp, KSPNormType normtype)
274: {
275:   PetscFunctionBegin;
278:   ksp->normtype = ksp->normtype_set = normtype;
279:   PetscFunctionReturn(PETSC_SUCCESS);
280: }

282: /*@
283:    KSPSetCheckNormIteration - Sets the first iteration at which the norm of the residual will be
284:      computed and used in the convergence test.

286:    Logically Collective

288:    Input Parameters:
289: +  ksp - Krylov solver context
290: -  it  - use -1 to check at all iterations

292:    Notes:
293:    Currently only works with `KSPCG`, `KSPBCGS` and `KSPIBCGS`

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

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

301: .seealso: `KSP`, `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSPConvergedSkip()`, `KSPSetNormType()`
302: @*/
303: PetscErrorCode KSPSetCheckNormIteration(KSP ksp, PetscInt it)
304: {
305:   PetscFunctionBegin;
308:   ksp->chknorm = it;
309:   PetscFunctionReturn(PETSC_SUCCESS);
310: }

312: /*@
313:    KSPSetLagNorm - Lags the residual norm calculation so that it is computed as part of the `MPI_Allreduce()` for
314:    computing the inner products for the next iteration.  This can reduce communication costs at the expense of doing
315:    one additional iteration.

317:    Logically Collective

319:    Input Parameters:
320: +  ksp - Krylov solver context
321: -  flg - `PETSC_TRUE` or `PETSC_FALSE`

323:    Options Database Key:
324: .  -ksp_lag_norm - lag the calculated residual norm

326:    Level: advanced

328:    Notes:
329:    Currently only works with `KSPIBCGS`.

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

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

335: .seealso: `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSPConvergedSkip()`, `KSPSetNormType()`, `KSPSetCheckNormIteration()`
336: @*/
337: PetscErrorCode KSPSetLagNorm(KSP ksp, PetscBool flg)
338: {
339:   PetscFunctionBegin;
342:   ksp->lagnorm = flg;
343:   PetscFunctionReturn(PETSC_SUCCESS);
344: }

346: /*@
347:    KSPSetSupportedNorm - Sets a norm and preconditioner side supported by a `KSP`

349:    Logically Collective

351:    Input Parameters:
352: +  ksp - Krylov method
353: .  normtype - supported norm type
354: .  pcside - preconditioner side that can be used with this norm
355: -  priority - positive integer preference for this combination; larger values have higher priority

357:    Level: developer

359:    Note:
360:    This function should be called from the implementation files `KSPCreate_XXX()` to declare
361:    which norms and preconditioner sides are supported. Users should not need to call this
362:    function.

364: .seealso: `KSP`, `KSPNormType`, `PCSide`, `KSPSetNormType()`, `KSPSetPCSide()`
365: @*/
366: PetscErrorCode KSPSetSupportedNorm(KSP ksp, KSPNormType normtype, PCSide pcside, PetscInt priority)
367: {
368:   PetscFunctionBegin;
370:   ksp->normsupporttable[normtype][pcside] = priority;
371:   PetscFunctionReturn(PETSC_SUCCESS);
372: }

374: PetscErrorCode KSPNormSupportTableReset_Private(KSP ksp)
375: {
376:   PetscFunctionBegin;
377:   PetscCall(PetscMemzero(ksp->normsupporttable, sizeof(ksp->normsupporttable)));
378:   ksp->pc_side  = ksp->pc_side_set;
379:   ksp->normtype = ksp->normtype_set;
380:   PetscFunctionReturn(PETSC_SUCCESS);
381: }

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

387:   PetscFunctionBegin;
388:   best = 0;
389:   for (i = 0; i < KSP_NORM_MAX; i++) {
390:     for (j = 0; j < PC_SIDE_MAX; j++) {
391:       if ((ksp->normtype == KSP_NORM_DEFAULT || ksp->normtype == i) && (ksp->pc_side == PC_SIDE_DEFAULT || ksp->pc_side == j) && (ksp->normsupporttable[i][j] > best)) {
392:         best  = ksp->normsupporttable[i][j];
393:         ibest = i;
394:         jbest = j;
395:       }
396:     }
397:   }
398:   if (best < 1 && errorifnotsupported) {
399:     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);
400:     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]);
401:     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]);
402:     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]);
403:   }
404:   if (normtype) *normtype = (KSPNormType)ibest;
405:   if (pcside) *pcside = (PCSide)jbest;
406:   PetscFunctionReturn(PETSC_SUCCESS);
407: }

409: /*@
410:    KSPGetNormType - Gets the norm that is used for convergence testing.

412:    Not Collective

414:    Input Parameter:
415: .  ksp - Krylov solver context

417:    Output Parameter:
418: .  normtype - norm that is used for convergence testing

420:    Level: advanced

422: .seealso: `KSPNormType`, `KSPSetNormType()`, `KSPConvergedSkip()`
423: @*/
424: PetscErrorCode KSPGetNormType(KSP ksp, KSPNormType *normtype)
425: {
426:   PetscFunctionBegin;
429:   PetscCall(KSPSetUpNorms_Private(ksp, PETSC_TRUE, &ksp->normtype, &ksp->pc_side));
430:   *normtype = ksp->normtype;
431:   PetscFunctionReturn(PETSC_SUCCESS);
432: }

434: #if defined(PETSC_HAVE_SAWS)
435: #include <petscviewersaws.h>
436: #endif

438: /*@
439:    KSPSetOperators - Sets the matrix associated with the linear system
440:    and a (possibly) different one from which the preconditioner will be built

442:    Collective

444:    Input Parameters:
445: +  ksp - the `KSP` context
446: .  Amat - the matrix that defines the linear system
447: -  Pmat - the matrix to be used in constructing the preconditioner, usually the same as Amat.

449:     Level: beginner

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

455:     All future calls to `KSPSetOperators()` must use the same size matrices!

457:     Passing a `NULL` for `Amat` or `Pmat` removes the matrix that is currently used.

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

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

472: .vb
473:          KSPGetOperators(ksp/pc,&mat,NULL); is equivalent to
474:            set size, type, etc of mat

476:          MatCreate(comm,&mat);
477:          KSP/PCSetOperators(ksp/pc,mat,mat);
478:          PetscObjectDereference((PetscObject)mat);
479:            set size, type, etc of mat

481:      and

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

486:          MatCreate(comm,&mat);
487:          MatCreate(comm,&pmat);
488:          KSP/PCSetOperators(ksp/pc,mat,pmat);
489:          PetscObjectDereference((PetscObject)mat);
490:          PetscObjectDereference((PetscObject)pmat);
491:            set size, type, etc of mat and pmat
492: .ve

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

503: .seealso: `KSP`, `Mat`, `KSPSolve()`, `KSPGetPC()`, `PCGetOperators()`, `PCSetOperators()`, `KSPGetOperators()`, `KSPSetComputeOperators()`, `KSPSetComputeInitialGuess()`, `KSPSetComputeRHS()`
504: @*/
505: PetscErrorCode KSPSetOperators(KSP ksp, Mat Amat, Mat Pmat)
506: {
507:   PetscFunctionBegin;
511:   if (Amat) PetscCheckSameComm(ksp, 1, Amat, 2);
512:   if (Pmat) PetscCheckSameComm(ksp, 1, Pmat, 3);
513:   if (!ksp->pc) PetscCall(KSPGetPC(ksp, &ksp->pc));
514:   PetscCall(PCSetOperators(ksp->pc, Amat, Pmat));
515:   if (ksp->setupstage == KSP_SETUP_NEWRHS) ksp->setupstage = KSP_SETUP_NEWMATRIX; /* so that next solve call will call PCSetUp() on new matrix */
516:   PetscFunctionReturn(PETSC_SUCCESS);
517: }

519: /*@
520:    KSPGetOperators - Gets the matrix associated with the linear system
521:    and a (possibly) different one used to construct the preconditioner.

523:    Collective

525:    Input Parameter:
526: .  ksp - the `KSP` context

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

532:     Level: intermediate

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

537: .seealso: `KSP`, `KSPSolve()`, `KSPGetPC()`, `PCGetOperators()`, `PCSetOperators()`, `KSPSetOperators()`, `KSPGetOperatorsSet()`
538: @*/
539: PetscErrorCode KSPGetOperators(KSP ksp, Mat *Amat, Mat *Pmat)
540: {
541:   PetscFunctionBegin;
543:   if (!ksp->pc) PetscCall(KSPGetPC(ksp, &ksp->pc));
544:   PetscCall(PCGetOperators(ksp->pc, Amat, Pmat));
545:   PetscFunctionReturn(PETSC_SUCCESS);
546: }

548: /*@C
549:    KSPGetOperatorsSet - Determines if the matrix associated with the linear system and
550:    possibly a different one associated with the preconditioner have been set in the `KSP`.

552:    Not Collective, though the results on all processes should be the same

554:    Input Parameter:
555: .  pc - the `KSP` context

557:    Output Parameters:
558: +  mat - the matrix associated with the linear system was set
559: -  pmat - matrix associated with the preconditioner was set, usually the same as `mat`

561:    Level: intermediate

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

567: .seealso: `KSP`, `PCSetOperators()`, `KSPGetOperators()`, `KSPSetOperators()`, `PCGetOperators()`, `PCGetOperatorsSet()`
568: @*/
569: PetscErrorCode KSPGetOperatorsSet(KSP ksp, PetscBool *mat, PetscBool *pmat)
570: {
571:   PetscFunctionBegin;
573:   if (!ksp->pc) PetscCall(KSPGetPC(ksp, &ksp->pc));
574:   PetscCall(PCGetOperatorsSet(ksp->pc, mat, pmat));
575:   PetscFunctionReturn(PETSC_SUCCESS);
576: }

578: /*@C
579:    KSPSetPreSolve - Sets a function that is called at the beginning of each `KSPSolve()`

581:    Logically Collective

583:    Input Parameters:
584: +   ksp - the solver object
585: .   presolve - the function to call before the solve
586: -   prectx - any context needed by the function

588:    Calling sequence of `presolve`:
589: $  PetscErrorCode func(KSP ksp, Vec rhs, Vec x, void *ctx)
590: +  ksp - the `KSP` context
591: .  rhs - the right-hand side vector
592: .  x - the solution vector
593: -  ctx - optional user-provided context

595:    Level: developer

597: .seealso: `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPSetPostSolve()`, `PCEISENSTAT`
598: @*/
599: PetscErrorCode KSPSetPreSolve(KSP ksp, PetscErrorCode (*presolve)(KSP, Vec, Vec, void *), void *prectx)
600: {
601:   PetscFunctionBegin;
603:   ksp->presolve = presolve;
604:   ksp->prectx   = prectx;
605:   PetscFunctionReturn(PETSC_SUCCESS);
606: }

608: /*@C
609:    KSPSetPostSolve - Sets a function that is called at the end of each `KSPSolve()` (whether it converges or not)

611:    Logically Collective

613:    Input Parameters:
614: +   ksp - the solver object
615: .   postsolve - the function to call after the solve
616: -   postctx - any context needed by the function

618:    Calling sequence of `postsolve`:
619: $  PetscErrorCode func(KSP ksp, Vec rhs, Vec x, void *ctx)
620: +  ksp - the `KSP` context
621: .  rhs - the right-hand side vector
622: .  x - the solution vector
623: -  ctx - optional user-provided context

625:    Level: developer

627: .seealso: `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPSetPreSolve()`, `PCEISENSTAT`
628: @*/
629: PetscErrorCode KSPSetPostSolve(KSP ksp, PetscErrorCode (*postsolve)(KSP, Vec, Vec, void *), void *postctx)
630: {
631:   PetscFunctionBegin;
633:   ksp->postsolve = postsolve;
634:   ksp->postctx   = postctx;
635:   PetscFunctionReturn(PETSC_SUCCESS);
636: }

638: /*@
639:    KSPCreate - Creates the `KSP` context.

641:    Collective

643:    Input Parameter:
644: .  comm - MPI communicator

646:    Output Parameter:
647: .  ksp - location to put the `KSP` context

649:    Level: beginner

651:    Note:
652:    The default `KSPType` is `KSPGMRES` with a restart of 30, using modified Gram-Schmidt orthogonalization.

654: .seealso: [](ch_ksp), `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPGMRES`, `KSPType`
655: @*/
656: PetscErrorCode KSPCreate(MPI_Comm comm, KSP *inksp)
657: {
658:   KSP   ksp;
659:   void *ctx;

661:   PetscFunctionBegin;
663:   *inksp = NULL;
664:   PetscCall(KSPInitializePackage());

666:   PetscCall(PetscHeaderCreate(ksp, KSP_CLASSID, "KSP", "Krylov Method", "KSP", comm, KSPDestroy, KSPView));

668:   ksp->max_it  = 10000;
669:   ksp->pc_side = ksp->pc_side_set = PC_SIDE_DEFAULT;
670:   ksp->rtol                       = 1.e-5;
671: #if defined(PETSC_USE_REAL_SINGLE)
672:   ksp->abstol = 1.e-25;
673: #else
674:   ksp->abstol = 1.e-50;
675: #endif
676:   ksp->divtol = 1.e4;

678:   ksp->chknorm  = -1;
679:   ksp->normtype = ksp->normtype_set = KSP_NORM_DEFAULT;
680:   ksp->rnorm                        = 0.0;
681:   ksp->its                          = 0;
682:   ksp->guess_zero                   = PETSC_TRUE;
683:   ksp->calc_sings                   = PETSC_FALSE;
684:   ksp->res_hist                     = NULL;
685:   ksp->res_hist_alloc               = NULL;
686:   ksp->res_hist_len                 = 0;
687:   ksp->res_hist_max                 = 0;
688:   ksp->res_hist_reset               = PETSC_TRUE;
689:   ksp->err_hist                     = NULL;
690:   ksp->err_hist_alloc               = NULL;
691:   ksp->err_hist_len                 = 0;
692:   ksp->err_hist_max                 = 0;
693:   ksp->err_hist_reset               = PETSC_TRUE;
694:   ksp->numbermonitors               = 0;
695:   ksp->numberreasonviews            = 0;
696:   ksp->setfromoptionscalled         = 0;
697:   ksp->nmax                         = PETSC_DECIDE;

699:   PetscCall(KSPConvergedDefaultCreate(&ctx));
700:   PetscCall(KSPSetConvergenceTest(ksp, KSPConvergedDefault, ctx, KSPConvergedDefaultDestroy));
701:   ksp->ops->buildsolution = KSPBuildSolutionDefault;
702:   ksp->ops->buildresidual = KSPBuildResidualDefault;

704:   ksp->vec_sol    = NULL;
705:   ksp->vec_rhs    = NULL;
706:   ksp->pc         = NULL;
707:   ksp->data       = NULL;
708:   ksp->nwork      = 0;
709:   ksp->work       = NULL;
710:   ksp->reason     = KSP_CONVERGED_ITERATING;
711:   ksp->setupstage = KSP_SETUP_NEW;

713:   PetscCall(KSPNormSupportTableReset_Private(ksp));

715:   *inksp = ksp;
716:   PetscFunctionReturn(PETSC_SUCCESS);
717: }

719: /*@C
720:    KSPSetType - Builds the `KSP` datastructure for a particular `KSPType`

722:    Logically Collective

724:    Input Parameters:
725: +  ksp  - the Krylov space context
726: -  type - a known method

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

731:    Level: intermediate

733:    Notes:
734:    See "petsc/include/petscksp.h" for available methods (for instance, `KSPCG` or `KSPGMRES`).

736:   Normally, it is best to use the `KSPSetFromOptions()` command and
737:   then set the `KSP` type from the options database rather than by using
738:   this routine.  Using the options database provides the user with
739:   maximum flexibility in evaluating the many different Krylov methods.
740:   The `KSPSetType()` routine is provided for those situations where it
741:   is necessary to set the iterative solver independently of the command
742:   line or options database.  This might be the case, for example, when
743:   the choice of iterative solver changes during the execution of the
744:   program, and the user's application is taking responsibility for
745:   choosing the appropriate method.  In other words, this routine is
746:   not for beginners.

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

751: .seealso: [](ch_ksp), `PCSetType()`, `KSPType`, `KSPRegister()`, `KSPCreate()`, `KSP`
752: @*/
753: PetscErrorCode KSPSetType(KSP ksp, KSPType type)
754: {
755:   PetscBool match;
756:   PetscErrorCode (*r)(KSP);

758:   PetscFunctionBegin;

762:   PetscCall(PetscObjectTypeCompare((PetscObject)ksp, type, &match));
763:   if (match) PetscFunctionReturn(PETSC_SUCCESS);

765:   PetscCall(PetscFunctionListFind(KSPList, type, &r));
766:   PetscCheck(r, PetscObjectComm((PetscObject)ksp), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unable to find requested KSP type %s", type);
767:   /* Destroy the previous private KSP context */
768:   PetscTryTypeMethod(ksp, destroy);
769:   ksp->ops->destroy = NULL;

771:   /* Reinitialize function pointers in KSPOps structure */
772:   PetscCall(PetscMemzero(ksp->ops, sizeof(struct _KSPOps)));
773:   ksp->ops->buildsolution = KSPBuildSolutionDefault;
774:   ksp->ops->buildresidual = KSPBuildResidualDefault;
775:   PetscCall(KSPNormSupportTableReset_Private(ksp));
776:   ksp->converged_neg_curve = PETSC_FALSE; // restore default
777:   ksp->setupnewmatrix      = PETSC_FALSE; // restore default (setup not called in case of new matrix)
778:   /* Call the KSPCreate_XXX routine for this particular Krylov solver */
779:   ksp->setupstage     = KSP_SETUP_NEW;
780:   ksp->guess_not_read = PETSC_FALSE; // restore default
781:   PetscCall((*r)(ksp));
782:   PetscCall(PetscObjectChangeTypeName((PetscObject)ksp, type));
783:   PetscFunctionReturn(PETSC_SUCCESS);
784: }

786: /*@C
787:    KSPGetType - Gets the `KSP` type as a string from the KSP object.

789:    Not Collective

791:    Input Parameter:
792: .  ksp - Krylov context

794:    Output Parameter:
795: .  name - name of the `KSP` method

797:    Level: intermediate

799: .seealso: [](ch_ksp), `KSPType`, `KSP`, `KSPSetType()`
800: @*/
801: PetscErrorCode KSPGetType(KSP ksp, KSPType *type)
802: {
803:   PetscFunctionBegin;
806:   *type = ((PetscObject)ksp)->type_name;
807:   PetscFunctionReturn(PETSC_SUCCESS);
808: }

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

813:    Not Collective

815:    Input Parameters:
816: +  sname - name of a new user-defined solver
817: -  function - routine to create method

819:    Level: advanced

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

824:    Sample usage:
825: .vb
826:    KSPRegister("my_solver", MySolverCreate);
827: .ve

829:    Then, your solver can be chosen with the procedural interface via
830: $    ` KSPSetType`(ksp, "my_solver")
831:    or at runtime via the option
832: $     -ksp_type my_solver

834: .seealso: [](ch_ksp), `KSP`, `KSPType`, `KSPSetType`, `KSPRegisterAll()`
835: @*/
836: PetscErrorCode KSPRegister(const char sname[], PetscErrorCode (*function)(KSP))
837: {
838:   PetscFunctionBegin;
839:   PetscCall(KSPInitializePackage());
840:   PetscCall(PetscFunctionListAdd(&KSPList, sname, function));
841:   PetscFunctionReturn(PETSC_SUCCESS);
842: }

844: PetscErrorCode KSPMonitorMakeKey_Internal(const char name[], PetscViewerType vtype, PetscViewerFormat format, char key[])
845: {
846:   PetscFunctionBegin;
847:   PetscCall(PetscStrncpy(key, name, PETSC_MAX_PATH_LEN));
848:   PetscCall(PetscStrlcat(key, ":", PETSC_MAX_PATH_LEN));
849:   PetscCall(PetscStrlcat(key, vtype, PETSC_MAX_PATH_LEN));
850:   PetscCall(PetscStrlcat(key, ":", PETSC_MAX_PATH_LEN));
851:   PetscCall(PetscStrlcat(key, PetscViewerFormats[format], PETSC_MAX_PATH_LEN));
852:   PetscFunctionReturn(PETSC_SUCCESS);
853: }

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

858:   Not Collective

860:   Input Parameters:
861: + name    - name of a new monitor routine
862: . vtype   - A `PetscViewerType` for the output
863: . format  - A `PetscViewerFormat` for the output
864: . monitor - Monitor routine
865: . create  - Creation routine, or `NULL`
866: - destroy - Destruction routine, or `NULL`

868:   Level: advanced

870:   Note:
871:   `KSPMonitorRegister()` may be called multiple times to add several user-defined monitors.

873:   Sample usage:
874: .vb
875:   KSPMonitorRegister("my_monitor", PETSCVIEWERASCII, PETSC_VIEWER_ASCII_INFO_DETAIL, MyMonitor, NULL, NULL);
876: .ve

878:   Then, your monitor can be chosen with the procedural interface via
879: $     KSPMonitorSetFromOptions(ksp, "-ksp_monitor_my_monitor", "my_monitor", NULL)
880:   or at runtime via the option
881: $     -ksp_monitor_my_monitor

883: .seealso: [](ch_ksp), `KSP`, `KSPMonitorSet()`, `KSPMonitorRegisterAll()`, `KSPMonitorSetFromOptions()`
884: @*/
885: PetscErrorCode KSPMonitorRegister(const char name[], PetscViewerType vtype, PetscViewerFormat format, PetscErrorCode (*monitor)(KSP, PetscInt, PetscReal, PetscViewerAndFormat *), PetscErrorCode (*create)(PetscViewer, PetscViewerFormat, void *, PetscViewerAndFormat **), PetscErrorCode (*destroy)(PetscViewerAndFormat **))
886: {
887:   char key[PETSC_MAX_PATH_LEN];

889:   PetscFunctionBegin;
890:   PetscCall(KSPInitializePackage());
891:   PetscCall(KSPMonitorMakeKey_Internal(name, vtype, format, key));
892:   PetscCall(PetscFunctionListAdd(&KSPMonitorList, key, monitor));
893:   if (create) PetscCall(PetscFunctionListAdd(&KSPMonitorCreateList, key, create));
894:   if (destroy) PetscCall(PetscFunctionListAdd(&KSPMonitorDestroyList, key, destroy));
895:   PetscFunctionReturn(PETSC_SUCCESS);
896: }