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;
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 binary with KSPView().
29: Collective on viewer
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: Notes:
39: The type is determined by the data in the file, any type set into the KSP before this call is ignored.
41: Notes for advanced users:
42: Most users should not need to know the details of the binary storage
43: format, since KSPLoad() and KSPView() completely hide these details.
44: But for anyone who's interested, the standard binary matrix storage
45: format is
46: .vb
47: has not yet been determined
48: .ve
50: .seealso: `PetscViewerBinaryOpen()`, `KSPView()`, `MatLoad()`, `VecLoad()`
51: @*/
52: PetscErrorCode KSPLoad(KSP newdm, PetscViewer viewer)
53: {
54: PetscBool isbinary;
55: PetscInt classid;
56: char type[256];
57: PC pc;
61: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary);
64: PetscViewerBinaryRead(viewer, &classid, 1, NULL, PETSC_INT);
66: PetscViewerBinaryRead(viewer, type, 256, NULL, PETSC_CHAR);
67: KSPSetType(newdm, type);
68: PetscTryTypeMethod(newdm, load, viewer);
69: KSPGetPC(newdm, &pc);
70: PCLoad(pc, viewer);
71: return 0;
72: }
74: #include <petscdraw.h>
75: #if defined(PETSC_HAVE_SAWS)
76: #include <petscviewersaws.h>
77: #endif
78: /*@C
79: KSPView - Prints the KSP data structure.
81: Collective on ksp
83: Input Parameters:
84: + ksp - the Krylov space context
85: - viewer - visualization context
87: Options Database Keys:
88: . -ksp_view - print the KSP data structure at the end of a KSPSolve call
90: Note:
91: The available visualization contexts include
92: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
93: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
94: output where only the first processor opens
95: the file. All other processors send their
96: data to the first processor to print.
98: The available formats include
99: + PETSC_VIEWER_DEFAULT - standard output (default)
100: - PETSC_VIEWER_ASCII_INFO_DETAIL - more verbose output for PCBJACOBI and PCASM
102: The user can open an alternative visualization context with
103: PetscViewerASCIIOpen() - output to a specified file.
105: In the debugger you can do "call KSPView(ksp,0)" to display the KSP. (The same holds for any PETSc object viewer).
107: Level: beginner
109: .seealso: `PCView()`, `PetscViewerASCIIOpen()`
110: @*/
111: PetscErrorCode KSPView(KSP ksp, PetscViewer viewer)
112: {
113: PetscBool iascii, isbinary, isdraw, isstring;
114: #if defined(PETSC_HAVE_SAWS)
115: PetscBool issaws;
116: #endif
119: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)ksp), &viewer);
123: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
124: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary);
125: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERDRAW, &isdraw);
126: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSTRING, &isstring);
127: #if defined(PETSC_HAVE_SAWS)
128: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSAWS, &issaws);
129: #endif
130: if (iascii) {
131: PetscObjectPrintClassNamePrefixType((PetscObject)ksp, viewer);
132: PetscViewerASCIIPushTab(viewer);
133: PetscTryTypeMethod(ksp, view, viewer);
134: PetscViewerASCIIPopTab(viewer);
135: if (ksp->guess_zero) {
136: PetscViewerASCIIPrintf(viewer, " maximum iterations=%" PetscInt_FMT ", initial guess is zero\n", ksp->max_it);
137: } else {
138: PetscViewerASCIIPrintf(viewer, " maximum iterations=%" PetscInt_FMT ", nonzero initial guess\n", ksp->max_it);
139: }
140: if (ksp->guess_knoll) PetscViewerASCIIPrintf(viewer, " using preconditioner applied to right hand side for initial guess\n");
141: PetscViewerASCIIPrintf(viewer, " tolerances: relative=%g, absolute=%g, divergence=%g\n", (double)ksp->rtol, (double)ksp->abstol, (double)ksp->divtol);
142: if (ksp->pc_side == PC_RIGHT) {
143: PetscViewerASCIIPrintf(viewer, " right preconditioning\n");
144: } else if (ksp->pc_side == PC_SYMMETRIC) {
145: PetscViewerASCIIPrintf(viewer, " symmetric preconditioning\n");
146: } else {
147: PetscViewerASCIIPrintf(viewer, " left preconditioning\n");
148: }
149: if (ksp->guess) {
150: PetscViewerASCIIPushTab(viewer);
151: KSPGuessView(ksp->guess, viewer);
152: PetscViewerASCIIPopTab(viewer);
153: }
154: if (ksp->dscale) PetscViewerASCIIPrintf(viewer, " diagonally scaled system\n");
155: PetscViewerASCIIPrintf(viewer, " using %s norm type for convergence test\n", KSPNormTypes[ksp->normtype]);
156: } else if (isbinary) {
157: PetscInt classid = KSP_FILE_CLASSID;
158: MPI_Comm comm;
159: PetscMPIInt rank;
160: char type[256];
162: PetscObjectGetComm((PetscObject)ksp, &comm);
163: MPI_Comm_rank(comm, &rank);
164: if (rank == 0) {
165: PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT);
166: PetscStrncpy(type, ((PetscObject)ksp)->type_name, 256);
167: PetscViewerBinaryWrite(viewer, type, 256, PETSC_CHAR);
168: }
169: PetscTryTypeMethod(ksp, view, viewer);
170: } else if (isstring) {
171: const char *type;
172: KSPGetType(ksp, &type);
173: PetscViewerStringSPrintf(viewer, " KSPType: %-7.7s", type);
174: PetscTryTypeMethod(ksp, view, viewer);
175: } else if (isdraw) {
176: PetscDraw draw;
177: char str[36];
178: PetscReal x, y, bottom, h;
179: PetscBool flg;
181: PetscViewerDrawGetDraw(viewer, 0, &draw);
182: PetscDrawGetCurrentPoint(draw, &x, &y);
183: PetscObjectTypeCompare((PetscObject)ksp, KSPPREONLY, &flg);
184: if (!flg) {
185: PetscStrncpy(str, "KSP: ", sizeof(str));
186: PetscStrlcat(str, ((PetscObject)ksp)->type_name, sizeof(str));
187: PetscDrawStringBoxed(draw, x, y, PETSC_DRAW_RED, PETSC_DRAW_BLACK, str, NULL, &h);
188: bottom = y - h;
189: } else {
190: bottom = y;
191: }
192: PetscDrawPushCurrentPoint(draw, x, bottom);
193: #if defined(PETSC_HAVE_SAWS)
194: } else if (issaws) {
195: PetscMPIInt rank;
196: const char *name;
198: PetscObjectGetName((PetscObject)ksp, &name);
199: MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
200: if (!((PetscObject)ksp)->amsmem && rank == 0) {
201: char dir[1024];
203: PetscObjectViewSAWs((PetscObject)ksp, viewer);
204: PetscSNPrintf(dir, 1024, "/PETSc/Objects/%s/its", name);
205: SAWs_Register, (dir, &ksp->its, 1, SAWs_READ, SAWs_INT);
206: if (!ksp->res_hist) KSPSetResidualHistory(ksp, NULL, PETSC_DECIDE, PETSC_TRUE);
207: PetscSNPrintf(dir, 1024, "/PETSc/Objects/%s/res_hist", name);
208: SAWs_Register, (dir, ksp->res_hist, 10, SAWs_READ, SAWs_DOUBLE);
209: }
210: #endif
211: } else PetscTryTypeMethod(ksp, view, viewer);
212: if (ksp->pc) PCView(ksp->pc, viewer);
213: if (isdraw) {
214: PetscDraw draw;
215: PetscViewerDrawGetDraw(viewer, 0, &draw);
216: PetscDrawPopCurrentPoint(draw);
217: }
218: return 0;
219: }
221: /*@C
222: KSPViewFromOptions - View from Options
224: Collective on KSP
226: Input Parameters:
227: + A - Krylov solver context
228: . obj - Optional object
229: - name - command line option
231: Level: intermediate
232: .seealso: `KSP`, `KSPView`, `PetscObjectViewFromOptions()`, `KSPCreate()`
233: @*/
234: PetscErrorCode KSPViewFromOptions(KSP A, PetscObject obj, const char name[])
235: {
237: PetscObjectViewFromOptions((PetscObject)A, obj, name);
238: return 0;
239: }
241: /*@
242: KSPSetNormType - Sets the norm that is used for convergence testing.
244: Logically Collective on ksp
246: Input Parameters:
247: + ksp - Krylov solver context
248: - normtype - one of
249: $ KSP_NORM_NONE - skips computing the norm, this should generally only be used if you are using
250: $ the Krylov method as a smoother with a fixed small number of iterations.
251: $ Implicitly sets KSPConvergedSkip() as KSP convergence test.
252: $ Note that certain algorithms such as KSPGMRES ALWAYS require the norm calculation,
253: $ for these methods the norms are still computed, they are just not used in
254: $ the convergence test.
255: $ KSP_NORM_PRECONDITIONED - the default for left preconditioned solves, uses the l2 norm
256: $ of the preconditioned residual P^{-1}(b - A x)
257: $ KSP_NORM_UNPRECONDITIONED - uses the l2 norm of the true b - Ax residual.
258: $ KSP_NORM_NATURAL - supported by KSPCG, KSPCR, KSPCGNE, KSPCGS
260: Options Database Key:
261: . -ksp_norm_type <none,preconditioned,unpreconditioned,natural> - set KSP norm type
263: Notes:
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 Notes:
269: Supported combinations of norm and preconditioner side are set using KSPSetSupportedNorm().
271: Level: advanced
273: .seealso: `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSPConvergedSkip()`, `KSPSetCheckNormIteration()`, `KSPSetPCSide()`, `KSPGetPCSide()`, `KSPNormType`
274: @*/
275: PetscErrorCode KSPSetNormType(KSP ksp, KSPNormType normtype)
276: {
279: ksp->normtype = ksp->normtype_set = normtype;
280: return 0;
281: }
283: /*@
284: KSPSetCheckNormIteration - Sets the first iteration at which the norm of the residual will be
285: computed and used in the convergence test.
287: Logically Collective on ksp
289: Input Parameters:
290: + ksp - Krylov solver context
291: - it - use -1 to check at all iterations
293: Notes:
294: Currently only works with KSPCG, KSPBCGS and KSPIBCGS
296: Use KSPSetNormType(ksp,KSP_NORM_NONE) to never check the norm
298: On steps where the norm is not computed, the previous norm is still in the variable, so if you run with, for example,
299: -ksp_monitor the residual norm will appear to be unchanged for several iterations (though it is not really unchanged).
300: Level: advanced
302: .seealso: `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSPConvergedSkip()`, `KSPSetNormType()`
303: @*/
304: PetscErrorCode KSPSetCheckNormIteration(KSP ksp, PetscInt it)
305: {
308: ksp->chknorm = it;
309: return 0;
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 on ksp
319: Input Parameters:
320: + ksp - Krylov solver context
321: - flg - PETSC_TRUE or PETSC_FALSE
323: Options Database Keys:
324: . -ksp_lag_norm - lag the calculated residual norm
326: Notes:
327: Currently only works with KSPIBCGS.
329: Use KSPSetNormType(ksp,KSP_NORM_NONE) to never check the norm
331: If you lag the norm and run with, for example, -ksp_monitor, the residual norm reported will be the lagged one.
332: Level: advanced
334: .seealso: `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSPConvergedSkip()`, `KSPSetNormType()`, `KSPSetCheckNormIteration()`
335: @*/
336: PetscErrorCode KSPSetLagNorm(KSP ksp, PetscBool flg)
337: {
340: ksp->lagnorm = flg;
341: return 0;
342: }
344: /*@
345: KSPSetSupportedNorm - Sets a norm and preconditioner side supported by a KSP
347: Logically Collective
349: Input Parameters:
350: + ksp - Krylov method
351: . normtype - supported norm type
352: . pcside - preconditioner side that can be used with this norm
353: - priority - positive integer preference for this combination; larger values have higher priority
355: Level: developer
357: Notes:
358: This function should be called from the implementation files KSPCreate_XXX() to declare
359: which norms and preconditioner sides are supported. Users should not need to call this
360: function.
362: .seealso: `KSPSetNormType()`, `KSPSetPCSide()`
363: @*/
364: PetscErrorCode KSPSetSupportedNorm(KSP ksp, KSPNormType normtype, PCSide pcside, PetscInt priority)
365: {
367: ksp->normsupporttable[normtype][pcside] = priority;
368: return 0;
369: }
371: PetscErrorCode KSPNormSupportTableReset_Private(KSP ksp)
372: {
373: PetscMemzero(ksp->normsupporttable, sizeof(ksp->normsupporttable));
374: ksp->pc_side = ksp->pc_side_set;
375: ksp->normtype = ksp->normtype_set;
376: return 0;
377: }
379: PetscErrorCode KSPSetUpNorms_Private(KSP ksp, PetscBool errorifnotsupported, KSPNormType *normtype, PCSide *pcside)
380: {
381: PetscInt i, j, best, ibest = 0, jbest = 0;
383: best = 0;
384: for (i = 0; i < KSP_NORM_MAX; i++) {
385: for (j = 0; j < PC_SIDE_MAX; j++) {
386: if ((ksp->normtype == KSP_NORM_DEFAULT || ksp->normtype == i) && (ksp->pc_side == PC_SIDE_DEFAULT || ksp->pc_side == j) && (ksp->normsupporttable[i][j] > best)) {
387: best = ksp->normsupporttable[i][j];
388: ibest = i;
389: jbest = j;
390: }
391: }
392: }
393: if (best < 1 && errorifnotsupported) {
397: 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]);
398: }
399: if (normtype) *normtype = (KSPNormType)ibest;
400: if (pcside) *pcside = (PCSide)jbest;
401: return 0;
402: }
404: /*@
405: KSPGetNormType - Gets the norm that is used for convergence testing.
407: Not Collective
409: Input Parameter:
410: . ksp - Krylov solver context
412: Output Parameter:
413: . normtype - norm that is used for convergence testing
415: Level: advanced
417: .seealso: `KSPNormType`, `KSPSetNormType()`, `KSPConvergedSkip()`
418: @*/
419: PetscErrorCode KSPGetNormType(KSP ksp, KSPNormType *normtype)
420: {
423: KSPSetUpNorms_Private(ksp, PETSC_TRUE, &ksp->normtype, &ksp->pc_side);
424: *normtype = ksp->normtype;
425: return 0;
426: }
428: #if defined(PETSC_HAVE_SAWS)
429: #include <petscviewersaws.h>
430: #endif
432: /*@
433: KSPSetOperators - Sets the matrix associated with the linear system
434: and a (possibly) different one associated with the preconditioner.
436: Collective on ksp
438: Input Parameters:
439: + ksp - the KSP context
440: . Amat - the matrix that defines the linear system
441: - Pmat - the matrix to be used in constructing the preconditioner, usually the same as Amat.
443: Notes:
445: If you know the operator Amat has a null space you can use MatSetNullSpace() and MatSetTransposeNullSpace() to supply the null
446: space to Amat and the KSP solvers will automatically use that null space as needed during the solution process.
448: All future calls to KSPSetOperators() must use the same size matrices!
450: Passing a NULL for Amat or Pmat removes the matrix that is currently used.
452: If you wish to replace either Amat or Pmat but leave the other one untouched then
453: first call KSPGetOperators() to get the one you wish to keep, call PetscObjectReference()
454: on it and then pass it back in in your call to KSPSetOperators().
456: Level: beginner
458: Alternative usage: If the operators have NOT been set with KSP/PCSetOperators() then the operators
459: are created in PC and returned to the user. In this case, if both operators
460: mat and pmat are requested, two DIFFERENT operators will be returned. If
461: only one is requested both operators in the PC will be the same (i.e. as
462: if one had called KSP/PCSetOperators() with the same argument for both Mats).
463: The user must set the sizes of the returned matrices and their type etc just
464: as if the user created them with MatCreate(). For example,
466: $ KSP/PCGetOperators(ksp/pc,&mat,NULL); is equivalent to
467: $ set size, type, etc of mat
469: $ MatCreate(comm,&mat);
470: $ KSP/PCSetOperators(ksp/pc,mat,mat);
471: $ PetscObjectDereference((PetscObject)mat);
472: $ set size, type, etc of mat
474: and
476: $ KSP/PCGetOperators(ksp/pc,&mat,&pmat); is equivalent to
477: $ set size, type, etc of mat and pmat
479: $ MatCreate(comm,&mat);
480: $ MatCreate(comm,&pmat);
481: $ KSP/PCSetOperators(ksp/pc,mat,pmat);
482: $ PetscObjectDereference((PetscObject)mat);
483: $ PetscObjectDereference((PetscObject)pmat);
484: $ set size, type, etc of mat and pmat
486: The rational for this support is so that when creating a TS, SNES, or KSP the hierarchy
487: of underlying objects (i.e. SNES, KSP, PC, Mat) and their livespans can be completely
488: managed by the top most level object (i.e. the TS, SNES, or KSP). Another way to look
489: at this is when you create a SNES you do not NEED to create a KSP and attach it to
490: the SNES object (the SNES object manages it for you). Similarly when you create a KSP
491: you do not need to attach a PC to it (the KSP object manages the PC object for you).
492: Thus, why should YOU have to create the Mat and attach it to the SNES/KSP/PC, when
493: it can be created for you?
495: .seealso: `KSPSolve()`, `KSPGetPC()`, `PCGetOperators()`, `PCSetOperators()`, `KSPGetOperators()`, `KSPSetComputeOperators()`, `KSPSetComputeInitialGuess()`, `KSPSetComputeRHS()`
496: @*/
497: PetscErrorCode KSPSetOperators(KSP ksp, Mat Amat, Mat Pmat)
498: {
504: if (!ksp->pc) KSPGetPC(ksp, &ksp->pc);
505: PCSetOperators(ksp->pc, Amat, Pmat);
506: if (ksp->setupstage == KSP_SETUP_NEWRHS) ksp->setupstage = KSP_SETUP_NEWMATRIX; /* so that next solve call will call PCSetUp() on new matrix */
507: return 0;
508: }
510: /*@
511: KSPGetOperators - Gets the matrix associated with the linear system
512: and a (possibly) different one associated with the preconditioner.
514: Collective on ksp
516: Input Parameter:
517: . ksp - the KSP context
519: Output Parameters:
520: + Amat - the matrix that defines the linear system
521: - Pmat - the matrix to be used in constructing the preconditioner, usually the same as Amat.
523: Level: intermediate
525: Notes:
526: DOES NOT increase the reference counts of the matrix, so you should NOT destroy them.
528: .seealso: `KSPSolve()`, `KSPGetPC()`, `PCGetOperators()`, `PCSetOperators()`, `KSPSetOperators()`, `KSPGetOperatorsSet()`
529: @*/
530: PetscErrorCode KSPGetOperators(KSP ksp, Mat *Amat, Mat *Pmat)
531: {
533: if (!ksp->pc) KSPGetPC(ksp, &ksp->pc);
534: PCGetOperators(ksp->pc, Amat, Pmat);
535: return 0;
536: }
538: /*@C
539: KSPGetOperatorsSet - Determines if the matrix associated with the linear system and
540: possibly a different one associated with the preconditioner have been set in the KSP.
542: Not collective, though the results on all processes should be the same
544: Input Parameter:
545: . pc - the KSP context
547: Output Parameters:
548: + mat - the matrix associated with the linear system was set
549: - pmat - matrix associated with the preconditioner was set, usually the same
551: Level: intermediate
553: .seealso: `PCSetOperators()`, `KSPGetOperators()`, `KSPSetOperators()`, `PCGetOperators()`, `PCGetOperatorsSet()`
554: @*/
555: PetscErrorCode KSPGetOperatorsSet(KSP ksp, PetscBool *mat, PetscBool *pmat)
556: {
558: if (!ksp->pc) KSPGetPC(ksp, &ksp->pc);
559: PCGetOperatorsSet(ksp->pc, mat, pmat);
560: return 0;
561: }
563: /*@C
564: KSPSetPreSolve - Sets a function that is called before every KSPSolve() is started
566: Logically Collective on ksp
568: Input Parameters:
569: + ksp - the solver object
570: . presolve - the function to call before the solve
571: - prectx - any context needed by the function
573: Calling sequence of presolve:
574: $ func(KSP ksp,Vec rhs,Vec x,void *ctx)
576: + ksp - the KSP context
577: . rhs - the right-hand side vector
578: . x - the solution vector
579: - ctx - optional user-provided context
581: Level: developer
583: .seealso: `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPSetPostSolve()`
584: @*/
585: PetscErrorCode KSPSetPreSolve(KSP ksp, PetscErrorCode (*presolve)(KSP, Vec, Vec, void *), void *prectx)
586: {
588: ksp->presolve = presolve;
589: ksp->prectx = prectx;
590: return 0;
591: }
593: /*@C
594: KSPSetPostSolve - Sets a function that is called after every KSPSolve() completes (whether it converges or not)
596: Logically Collective on ksp
598: Input Parameters:
599: + ksp - the solver object
600: . postsolve - the function to call after the solve
601: - postctx - any context needed by the function
603: Level: developer
605: Calling sequence of postsolve:
606: $ func(KSP ksp,Vec rhs,Vec x,void *ctx)
608: + ksp - the KSP context
609: . rhs - the right-hand side vector
610: . x - the solution vector
611: - ctx - optional user-provided context
613: .seealso: `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPSetPreSolve()`
614: @*/
615: PetscErrorCode KSPSetPostSolve(KSP ksp, PetscErrorCode (*postsolve)(KSP, Vec, Vec, void *), void *postctx)
616: {
618: ksp->postsolve = postsolve;
619: ksp->postctx = postctx;
620: return 0;
621: }
623: /*@
624: KSPCreate - Creates the default KSP context.
626: Collective
628: Input Parameter:
629: . comm - MPI communicator
631: Output Parameter:
632: . ksp - location to put the KSP context
634: Notes:
635: The default KSP type is GMRES with a restart of 30, using modified Gram-Schmidt
636: orthogonalization.
638: Level: beginner
640: .seealso: `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`
641: @*/
642: PetscErrorCode KSPCreate(MPI_Comm comm, KSP *inksp)
643: {
644: KSP ksp;
645: void *ctx;
648: *inksp = NULL;
649: KSPInitializePackage();
651: PetscHeaderCreate(ksp, KSP_CLASSID, "KSP", "Krylov Method", "KSP", comm, KSPDestroy, KSPView);
653: ksp->max_it = 10000;
654: ksp->pc_side = ksp->pc_side_set = PC_SIDE_DEFAULT;
655: ksp->rtol = 1.e-5;
656: #if defined(PETSC_USE_REAL_SINGLE)
657: ksp->abstol = 1.e-25;
658: #else
659: ksp->abstol = 1.e-50;
660: #endif
661: ksp->divtol = 1.e4;
663: ksp->chknorm = -1;
664: ksp->normtype = ksp->normtype_set = KSP_NORM_DEFAULT;
665: ksp->rnorm = 0.0;
666: ksp->its = 0;
667: ksp->guess_zero = PETSC_TRUE;
668: ksp->calc_sings = PETSC_FALSE;
669: ksp->res_hist = NULL;
670: ksp->res_hist_alloc = NULL;
671: ksp->res_hist_len = 0;
672: ksp->res_hist_max = 0;
673: ksp->res_hist_reset = PETSC_TRUE;
674: ksp->err_hist = NULL;
675: ksp->err_hist_alloc = NULL;
676: ksp->err_hist_len = 0;
677: ksp->err_hist_max = 0;
678: ksp->err_hist_reset = PETSC_TRUE;
679: ksp->numbermonitors = 0;
680: ksp->numberreasonviews = 0;
681: ksp->setfromoptionscalled = 0;
682: ksp->nmax = PETSC_DECIDE;
684: KSPConvergedDefaultCreate(&ctx);
685: KSPSetConvergenceTest(ksp, KSPConvergedDefault, ctx, KSPConvergedDefaultDestroy);
686: ksp->ops->buildsolution = KSPBuildSolutionDefault;
687: ksp->ops->buildresidual = KSPBuildResidualDefault;
689: ksp->vec_sol = NULL;
690: ksp->vec_rhs = NULL;
691: ksp->pc = NULL;
692: ksp->data = NULL;
693: ksp->nwork = 0;
694: ksp->work = NULL;
695: ksp->reason = KSP_CONVERGED_ITERATING;
696: ksp->setupstage = KSP_SETUP_NEW;
698: KSPNormSupportTableReset_Private(ksp);
700: *inksp = ksp;
701: return 0;
702: }
704: /*@C
705: KSPSetType - Builds KSP for a particular solver.
707: Logically Collective on ksp
709: Input Parameters:
710: + ksp - the Krylov space context
711: - type - a known method
713: Options Database Key:
714: . -ksp_type <method> - Sets the method; use -help for a list
715: of available methods (for instance, cg or gmres)
717: Notes:
718: See "petsc/include/petscksp.h" for available methods (for instance,
719: KSPCG or KSPGMRES).
721: Normally, it is best to use the KSPSetFromOptions() command and
722: then set the KSP type from the options database rather than by using
723: this routine. Using the options database provides the user with
724: maximum flexibility in evaluating the many different Krylov methods.
725: The KSPSetType() routine is provided for those situations where it
726: is necessary to set the iterative solver independently of the command
727: line or options database. This might be the case, for example, when
728: the choice of iterative solver changes during the execution of the
729: program, and the user's application is taking responsibility for
730: choosing the appropriate method. In other words, this routine is
731: not for beginners.
733: Level: intermediate
735: Developer Note: KSPRegister() is used to add Krylov types to KSPList from which they
736: are accessed by KSPSetType().
738: .seealso: `PCSetType()`, `KSPType`, `KSPRegister()`, `KSPCreate()`
740: @*/
741: PetscErrorCode KSPSetType(KSP ksp, KSPType type)
742: {
743: PetscBool match;
744: PetscErrorCode (*r)(KSP);
749: PetscObjectTypeCompare((PetscObject)ksp, type, &match);
750: if (match) return 0;
752: PetscFunctionListFind(KSPList, type, &r);
754: /* Destroy the previous private KSP context */
755: PetscTryTypeMethod(ksp, destroy);
756: ksp->ops->destroy = NULL;
758: /* Reinitialize function pointers in KSPOps structure */
759: PetscMemzero(ksp->ops, sizeof(struct _KSPOps));
760: ksp->ops->buildsolution = KSPBuildSolutionDefault;
761: ksp->ops->buildresidual = KSPBuildResidualDefault;
762: KSPNormSupportTableReset_Private(ksp);
763: ksp->setupnewmatrix = PETSC_FALSE; // restore default (setup not called in case of new matrix)
764: /* Call the KSPCreate_XXX routine for this particular Krylov solver */
765: ksp->setupstage = KSP_SETUP_NEW;
766: (*r)(ksp);
767: PetscObjectChangeTypeName((PetscObject)ksp, type);
768: return 0;
769: }
771: /*@C
772: KSPGetType - Gets the KSP type as a string from the KSP object.
774: Not Collective
776: Input Parameter:
777: . ksp - Krylov context
779: Output Parameter:
780: . name - name of KSP method
782: Level: intermediate
784: .seealso: `KSPSetType()`
785: @*/
786: PetscErrorCode KSPGetType(KSP ksp, KSPType *type)
787: {
790: *type = ((PetscObject)ksp)->type_name;
791: return 0;
792: }
794: /*@C
795: KSPRegister - Adds a method to the Krylov subspace solver package.
797: Not Collective
799: Input Parameters:
800: + name_solver - name of a new user-defined solver
801: - routine_create - routine to create method context
803: Notes:
804: KSPRegister() may be called multiple times to add several user-defined solvers.
806: Sample usage:
807: .vb
808: KSPRegister("my_solver",MySolverCreate);
809: .ve
811: Then, your solver can be chosen with the procedural interface via
812: $ KSPSetType(ksp,"my_solver")
813: or at runtime via the option
814: $ -ksp_type my_solver
816: Level: advanced
818: .seealso: `KSPRegisterAll()`
819: @*/
820: PetscErrorCode KSPRegister(const char sname[], PetscErrorCode (*function)(KSP))
821: {
822: KSPInitializePackage();
823: PetscFunctionListAdd(&KSPList, sname, function);
824: return 0;
825: }
827: PetscErrorCode KSPMonitorMakeKey_Internal(const char name[], PetscViewerType vtype, PetscViewerFormat format, char key[])
828: {
829: PetscStrncpy(key, name, PETSC_MAX_PATH_LEN);
830: PetscStrlcat(key, ":", PETSC_MAX_PATH_LEN);
831: PetscStrlcat(key, vtype, PETSC_MAX_PATH_LEN);
832: PetscStrlcat(key, ":", PETSC_MAX_PATH_LEN);
833: PetscStrlcat(key, PetscViewerFormats[format], PETSC_MAX_PATH_LEN);
834: return 0;
835: }
837: /*@C
838: KSPMonitorRegister - Adds Krylov subspace solver monitor routine.
840: Not Collective
842: Input Parameters:
843: + name - name of a new monitor routine
844: . vtype - A PetscViewerType for the output
845: . format - A PetscViewerFormat for the output
846: . monitor - Monitor routine
847: . create - Creation routine, or NULL
848: - destroy - Destruction routine, or NULL
850: Notes:
851: KSPMonitorRegister() may be called multiple times to add several user-defined monitors.
853: Sample usage:
854: .vb
855: KSPMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
856: .ve
858: Then, your monitor can be chosen with the procedural interface via
859: $ KSPMonitorSetFromOptions(ksp,"-ksp_monitor_my_monitor","my_monitor",NULL)
860: or at runtime via the option
861: $ -ksp_monitor_my_monitor
863: Level: advanced
865: .seealso: `KSPMonitorRegisterAll()`
866: @*/
867: PetscErrorCode KSPMonitorRegister(const char name[], PetscViewerType vtype, PetscViewerFormat format, PetscErrorCode (*monitor)(KSP, PetscInt, PetscReal, PetscViewerAndFormat *), PetscErrorCode (*create)(PetscViewer, PetscViewerFormat, void *, PetscViewerAndFormat **), PetscErrorCode (*destroy)(PetscViewerAndFormat **))
868: {
869: char key[PETSC_MAX_PATH_LEN];
871: KSPInitializePackage();
872: KSPMonitorMakeKey_Internal(name, vtype, format, key);
873: PetscFunctionListAdd(&KSPMonitorList, key, monitor);
874: if (create) PetscFunctionListAdd(&KSPMonitorCreateList, key, create);
875: if (destroy) PetscFunctionListAdd(&KSPMonitorDestroyList, key, destroy);
876: return 0;
877: }