Actual source code: itcl.c
2: /*
3: Code for setting KSP options from the options database.
4: */
6: #include <petsc/private/kspimpl.h>
7: #include <petscdraw.h>
9: /*@C
10: KSPSetOptionsPrefix - Sets the prefix used for searching for all
11: `KSP` options in the database.
13: Logically Collective
15: Input Parameters:
16: + ksp - the Krylov context
17: - prefix - the prefix string to prepend to all `KSP` option requests
19: Level: advanced
21: Notes:
22: A hyphen (-) must NOT be given at the beginning of the prefix name.
23: The first character of all runtime options is AUTOMATICALLY the
24: hyphen.
26: For example, to distinguish between the runtime options for two
27: different `KSP` contexts, one could call
28: .vb
29: KSPSetOptionsPrefix(ksp1,"sys1_")
30: KSPSetOptionsPrefix(ksp2,"sys2_")
31: .ve
33: This would enable use of different options for each system, such as
34: .vb
35: -sys1_ksp_type gmres -sys1_ksp_rtol 1.e-3
36: -sys2_ksp_type bcgs -sys2_ksp_rtol 1.e-4
37: .ve
39: .seealso: [](ch_ksp), `KSP`, `KSPAppendOptionsPrefix()`, `KSPGetOptionsPrefix()`, `KSPSetFromOptions()`
40: @*/
41: PetscErrorCode KSPSetOptionsPrefix(KSP ksp, const char prefix[])
42: {
43: PetscFunctionBegin;
45: if (!ksp->pc) PetscCall(KSPGetPC(ksp, &ksp->pc));
46: PetscCall(PCSetOptionsPrefix(ksp->pc, prefix));
47: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)ksp, prefix));
48: PetscFunctionReturn(PETSC_SUCCESS);
49: }
51: /*@C
52: KSPAppendOptionsPrefix - Appends to the prefix used for searching for all
53: `KSP` options in the database.
55: Logically Collective
57: Input Parameters:
58: + ksp - the Krylov context
59: - prefix - the prefix string to prepend to all `KSP` option requests
61: Notes:
62: A hyphen (-) must NOT be given at the beginning of the prefix name.
63: The first character of all runtime options is AUTOMATICALLY the hyphen.
65: Level: advanced
67: .seealso: [](ch_ksp), `KSP`, `KSPSetOptionsPrefix()`, `KSPGetOptionsPrefix()`, `KSPSetFromOptions()`
68: @*/
69: PetscErrorCode KSPAppendOptionsPrefix(KSP ksp, const char prefix[])
70: {
71: PetscFunctionBegin;
73: if (!ksp->pc) PetscCall(KSPGetPC(ksp, &ksp->pc));
74: PetscCall(PCAppendOptionsPrefix(ksp->pc, prefix));
75: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)ksp, prefix));
76: PetscFunctionReturn(PETSC_SUCCESS);
77: }
79: /*@
80: KSPSetUseFischerGuess - Use the Paul Fischer algorithm or its variants to compute initial guesses for a set of solves with related right hand sides
82: Logically Collective
84: Input Parameters:
85: + ksp - the Krylov context
86: . model - use model 1, model 2, model 3, or any other number to turn it off
87: - size - size of subspace used to generate initial guess
89: Options Database Key:
90: . -ksp_fischer_guess <model,size> - uses the Fischer initial guess generator for repeated linear solves
92: Level: advanced
94: .seealso: [](ch_ksp), `KSP`, `KSPSetOptionsPrefix()`, `KSPAppendOptionsPrefix()`, `KSPSetUseFischerGuess()`, `KSPSetGuess()`, `KSPGetGuess()`, `KSPGuess`
95: @*/
96: PetscErrorCode KSPSetUseFischerGuess(KSP ksp, PetscInt model, PetscInt size)
97: {
98: KSPGuess guess;
100: PetscFunctionBegin;
104: PetscCall(KSPGetGuess(ksp, &guess));
105: PetscCall(KSPGuessSetType(guess, KSPGUESSFISCHER));
106: PetscCall(KSPGuessFischerSetModel(guess, model, size));
107: PetscFunctionReturn(PETSC_SUCCESS);
108: }
110: /*@
111: KSPSetGuess - Set the initial guess object
113: Logically Collective
115: Input Parameters:
116: + ksp - the Krylov context
117: - guess - the object created with `KSPGuessCreate()`
119: Level: advanced
121: Notes:
122: this allows a single `KSP` to be used with several different initial guess generators (likely for different linear
123: solvers, see `KSPSetPC()`).
125: This increases the reference count of the guess object, you must destroy the object with `KSPGuessDestroy()`
126: before the end of the program.
128: .seealso: [](ch_ksp), `KSP`, `KSPGuess`, `KSPSetOptionsPrefix()`, `KSPAppendOptionsPrefix()`, `KSPSetUseFischerGuess()`, `KSPGetGuess()`
129: @*/
130: PetscErrorCode KSPSetGuess(KSP ksp, KSPGuess guess)
131: {
132: PetscFunctionBegin;
135: PetscCall(PetscObjectReference((PetscObject)guess));
136: PetscCall(KSPGuessDestroy(&ksp->guess));
137: ksp->guess = guess;
138: ksp->guess->ksp = ksp;
139: PetscFunctionReturn(PETSC_SUCCESS);
140: }
142: /*@
143: KSPGetGuess - Gets the initial guess generator for the `KSP`.
145: Not Collective
147: Input Parameter:
148: . ksp - the Krylov context
150: Output Parameter:
151: . guess - the object
153: Level: developer
155: .seealso: [](ch_ksp), `KSPGuess`, `KSP`, `KSPSetOptionsPrefix()`, `KSPAppendOptionsPrefix()`, `KSPSetUseFischerGuess()`, `KSPSetGuess()`
156: @*/
157: PetscErrorCode KSPGetGuess(KSP ksp, KSPGuess *guess)
158: {
159: PetscFunctionBegin;
162: if (!ksp->guess) {
163: const char *prefix;
165: PetscCall(KSPGuessCreate(PetscObjectComm((PetscObject)ksp), &ksp->guess));
166: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix));
167: if (prefix) PetscCall(PetscObjectSetOptionsPrefix((PetscObject)ksp->guess, prefix));
168: ksp->guess->ksp = ksp;
169: }
170: *guess = ksp->guess;
171: PetscFunctionReturn(PETSC_SUCCESS);
172: }
174: /*@C
175: KSPGetOptionsPrefix - Gets the prefix used for searching for all
176: `KSP` options in the database.
178: Not Collective
180: Input Parameter:
181: . ksp - the Krylov context
183: Output Parameter:
184: . prefix - pointer to the prefix string used is returned
186: Notes:
187: On the fortran side, the user should pass in a string 'prefix' of
188: sufficient length to hold the prefix.
190: Level: advanced
192: .seealso: [](ch_ksp), `KSP`, `KSPSetFromOptions()`, `KSPSetOptionsPrefix()`, `KSPAppendOptionsPrefix()`
193: @*/
194: PetscErrorCode KSPGetOptionsPrefix(KSP ksp, const char *prefix[])
195: {
196: PetscFunctionBegin;
198: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)ksp, prefix));
199: PetscFunctionReturn(PETSC_SUCCESS);
200: }
202: static PetscErrorCode PetscViewerAndFormatCreate_Internal(PetscViewer viewer, PetscViewerFormat format, void *ctx, PetscViewerAndFormat **vf)
203: {
204: PetscFunctionBegin;
205: PetscCall(PetscViewerAndFormatCreate(viewer, format, vf));
206: (*vf)->data = ctx;
207: PetscFunctionReturn(PETSC_SUCCESS);
208: }
210: /*@C
211: KSPMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user in the options database
213: Collective
215: Input Parameters:
216: + ksp - `KSP` object you wish to monitor
217: . opt - the command line option for this monitor
218: . name - the monitor type one is seeking
219: - ctx - An optional user context for the monitor, or NULL
221: Level: developer
223: .seealso: [](ch_ksp), `KSPMonitorRegister()`, `KSPMonitorSet()`, `PetscOptionsGetViewer()`, `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`,
224: `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`
225: `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`, `PetscOptionsBool()`,
226: `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`,
227: `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`,
228: `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`,
229: `PetscOptionsFList()`, `PetscOptionsEList()`
230: @*/
231: PetscErrorCode KSPMonitorSetFromOptions(KSP ksp, const char opt[], const char name[], void *ctx)
232: {
233: PetscErrorCode (*mfunc)(KSP, PetscInt, PetscReal, void *);
234: PetscErrorCode (*cfunc)(PetscViewer, PetscViewerFormat, void *, PetscViewerAndFormat **);
235: PetscErrorCode (*dfunc)(PetscViewerAndFormat **);
236: PetscViewerAndFormat *vf;
237: PetscViewer viewer;
238: PetscViewerFormat format;
239: PetscViewerType vtype;
240: char key[PETSC_MAX_PATH_LEN];
241: PetscBool all, flg;
242: const char *prefix = NULL;
244: PetscFunctionBegin;
245: PetscCall(PetscStrcmp(opt, "-all_ksp_monitor", &all));
246: if (!all) PetscCall(PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix));
247: PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)ksp), ((PetscObject)ksp)->options, prefix, opt, &viewer, &format, &flg));
248: if (!flg) PetscFunctionReturn(PETSC_SUCCESS);
250: PetscCall(PetscViewerGetType(viewer, &vtype));
251: PetscCall(KSPMonitorMakeKey_Internal(name, vtype, format, key));
252: PetscCall(PetscFunctionListFind(KSPMonitorList, key, &mfunc));
253: PetscCall(PetscFunctionListFind(KSPMonitorCreateList, key, &cfunc));
254: PetscCall(PetscFunctionListFind(KSPMonitorDestroyList, key, &dfunc));
255: if (!cfunc) cfunc = PetscViewerAndFormatCreate_Internal;
256: if (!dfunc) dfunc = PetscViewerAndFormatDestroy;
258: PetscCall((*cfunc)(viewer, format, ctx, &vf));
259: PetscCall(PetscObjectDereference((PetscObject)viewer));
260: PetscCall(KSPMonitorSet(ksp, mfunc, vf, (PetscErrorCode(*)(void **))dfunc));
261: PetscFunctionReturn(PETSC_SUCCESS);
262: }
264: /*@
265: KSPSetFromOptions - Sets `KSP` options from the options database.
266: This routine must be called before `KSPSetUp()` if the user is to be
267: allowed to set the Krylov type.
269: Collective
271: Input Parameter:
272: . ksp - the Krylov space context
274: Options Database Keys:
275: + -ksp_max_it - maximum number of linear iterations
276: . -ksp_rtol rtol - relative tolerance used in default determination of convergence, i.e.
277: if residual norm decreases by this factor than convergence is declared
278: . -ksp_atol abstol - absolute tolerance used in default convergence test, i.e. if residual
279: norm is less than this then convergence is declared
280: . -ksp_divtol tol - if residual norm increases by this factor than divergence is declared
281: . -ksp_converged_use_initial_residual_norm - see `KSPConvergedDefaultSetUIRNorm()`
282: . -ksp_converged_use_min_initial_residual_norm - see `KSPConvergedDefaultSetUMIRNorm()`
283: . -ksp_converged_maxits - see `KSPConvergedDefaultSetConvergedMaxits()`
284: . -ksp_norm_type - none - skip norms used in convergence tests (useful only when not using
285: convergence test (say you always want to run with 5 iterations) to
286: save on communication overhead
287: preconditioned - default for left preconditioning
288: unpreconditioned - see `KSPSetNormType()`
289: natural - see `KSPSetNormType()`
290: . -ksp_check_norm_iteration it - do not compute residual norm until iteration number it (does compute at 0th iteration)
291: works only for `KSPBCGS`, `KSPIBCGS` and and `KSPCG`
292: . -ksp_lag_norm - compute the norm of the residual for the ith iteration on the i+1 iteration; this means that one can use
293: the norm of the residual for convergence test WITHOUT an extra `MPI_Allreduce()` limiting global synchronizations.
294: This will require 1 more iteration of the solver than usual.
295: . -ksp_guess_type - Type of initial guess generator for repeated linear solves
296: . -ksp_fischer_guess <model,size> - uses the Fischer initial guess generator for repeated linear solves
297: . -ksp_constant_null_space - assume the operator (matrix) has the constant vector in its null space
298: . -ksp_test_null_space - tests the null space set with `MatSetNullSpace()` to see if it truly is a null space
299: . -ksp_knoll - compute initial guess by applying the preconditioner to the right hand side
300: . -ksp_monitor_cancel - cancel all previous convergene monitor routines set
301: . -ksp_monitor - print residual norm at each iteration
302: . -ksp_monitor draw::draw_lg - plot residual norm at each iteration
303: . -ksp_monitor_true_residual - print true residual norm at each iteration
304: . -all_ksp_monitor <optional filename> - print residual norm at each iteration for ALL KSP solves, regardless of their prefix. This is
305: useful for `PCFIELDSPLIT`, `PCMG`, etc that have inner solvers and you wish to track the convergence of all the solvers
306: . -ksp_monitor_solution [ascii binary or draw][:filename][:format option] - plot solution at each iteration
307: . -ksp_monitor_singular_value - monitor extreme singular values at each iteration
308: . -ksp_converged_reason - view the convergence state at the end of the solve
309: . -ksp_use_explicittranspose - transpose the system explicitly in KSPSolveTranspose
310: . -ksp_error_if_not_converged - stop the program as soon as an error is detected in a `KSPSolve()`, `KSP_DIVERGED_ITS` is not treated as an error on inner solves
311: - -ksp_converged_rate - view the convergence rate at the end of the solve
313: Notes:
314: To see all options, run your program with the -help option
315: or consult [](ch_ksp)
317: Level: beginner
319: .seealso: [](ch_ksp), `KSP`, `KSPSetOptionsPrefix()`, `KSPResetFromOptions()`, `KSPSetUseFischerGuess()`
320: @*/
321: PetscErrorCode KSPSetFromOptions(KSP ksp)
322: {
323: const char *convtests[] = {"default", "skip", "lsqr"}, *prefix;
324: char type[256], guesstype[256], monfilename[PETSC_MAX_PATH_LEN];
325: PetscBool flg, flag, reuse, set;
326: PetscInt indx, model[2] = {0, 0}, nmax;
327: KSPNormType normtype;
328: PCSide pcside;
329: void *ctx;
330: MPI_Comm comm;
332: PetscFunctionBegin;
334: PetscCall(PetscObjectGetComm((PetscObject)ksp, &comm));
335: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)ksp, &prefix));
336: if (!ksp->skippcsetfromoptions) {
337: if (!ksp->pc) PetscCall(KSPGetPC(ksp, &ksp->pc));
338: PetscCall(PCSetFromOptions(ksp->pc));
339: }
341: PetscCall(KSPRegisterAll());
342: PetscObjectOptionsBegin((PetscObject)ksp);
343: PetscCall(PetscOptionsFList("-ksp_type", "Krylov method", "KSPSetType", KSPList, (char *)(((PetscObject)ksp)->type_name ? ((PetscObject)ksp)->type_name : KSPGMRES), type, 256, &flg));
344: if (flg) PetscCall(KSPSetType(ksp, type));
345: /*
346: Set the type if it was never set.
347: */
348: if (!((PetscObject)ksp)->type_name) PetscCall(KSPSetType(ksp, KSPGMRES));
350: PetscCall(KSPResetViewers(ksp));
352: /* Cancels all monitors hardwired into code before call to KSPSetFromOptions() */
353: PetscCall(PetscOptionsBool("-ksp_monitor_cancel", "Remove any hardwired monitor routines", "KSPMonitorCancel", PETSC_FALSE, &flg, &set));
354: if (set && flg) PetscCall(KSPMonitorCancel(ksp));
355: PetscCall(KSPMonitorSetFromOptions(ksp, "-ksp_monitor", "preconditioned_residual", NULL));
356: PetscCall(KSPMonitorSetFromOptions(ksp, "-ksp_monitor_short", "preconditioned_residual_short", NULL));
357: PetscCall(KSPMonitorSetFromOptions(ksp, "-all_ksp_monitor", "preconditioned_residual", NULL));
358: PetscCall(KSPMonitorSetFromOptions(ksp, "-ksp_monitor_range", "preconditioned_residual_range", NULL));
359: PetscCall(KSPMonitorSetFromOptions(ksp, "-ksp_monitor_true_residual", "true_residual", NULL));
360: PetscCall(KSPMonitorSetFromOptions(ksp, "-ksp_monitor_max", "true_residual_max", NULL));
361: PetscCall(KSPMonitorSetFromOptions(ksp, "-ksp_monitor_solution", "solution", NULL));
362: PetscCall(KSPMonitorSetFromOptions(ksp, "-ksp_monitor_singular_value", "singular_value", ksp));
363: PetscCall(KSPMonitorSetFromOptions(ksp, "-ksp_monitor_error", "error", ksp));
364: PetscCall(PetscOptionsBool("-ksp_monitor_pause_final", "Pauses all draw monitors at the final iterate", "KSPMonitorPauseFinal_Internal", PETSC_FALSE, &ksp->pauseFinal, NULL));
365: PetscCall(PetscOptionsBool("-ksp_initial_guess_nonzero", "Use the contents of the solution vector for initial guess", "KSPSetInitialNonzero", ksp->guess_zero ? PETSC_FALSE : PETSC_TRUE, &flag, &flg));
366: if (flg) PetscCall(KSPSetInitialGuessNonzero(ksp, flag));
368: PetscCall(PetscObjectTypeCompare((PetscObject)ksp, KSPPREONLY, &flg));
369: if (flg) {
370: PetscCall(KSPGetReusePreconditioner(ksp, &reuse));
371: PetscCall(PetscOptionsBool("-ksp_reuse_preconditioner", "Use initial preconditioner and don't ever compute a new one", "KSPReusePreconditioner", reuse, &reuse, NULL));
372: PetscCall(KSPSetReusePreconditioner(ksp, reuse));
373: PetscCall(PetscOptionsBool("-ksp_error_if_not_converged", "Generate error if solver does not converge", "KSPSetErrorIfNotConverged", ksp->errorifnotconverged, &ksp->errorifnotconverged, NULL));
374: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view", &ksp->viewer, &ksp->format, &ksp->view));
375: flg = PETSC_FALSE;
376: PetscCall(PetscOptionsBool("-ksp_converged_reason_view_cancel", "Cancel all the converged reason view functions set using KSPConvergedReasonViewSet", "KSPConvergedReasonViewCancel", PETSC_FALSE, &flg, &set));
377: if (set && flg) PetscCall(KSPConvergedReasonViewCancel(ksp));
378: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_mat", &ksp->viewerMat, &ksp->formatMat, &ksp->viewMat));
379: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_pmat", &ksp->viewerPMat, &ksp->formatPMat, &ksp->viewPMat));
380: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_rhs", &ksp->viewerRhs, &ksp->formatRhs, &ksp->viewRhs));
381: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_solution", &ksp->viewerSol, &ksp->formatSol, &ksp->viewSol));
382: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_mat_explicit", &ksp->viewerMatExp, &ksp->formatMatExp, &ksp->viewMatExp));
383: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_final_residual", &ksp->viewerFinalRes, &ksp->formatFinalRes, &ksp->viewFinalRes));
384: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_preconditioned_operator_explicit", &ksp->viewerPOpExp, &ksp->formatPOpExp, &ksp->viewPOpExp));
385: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_diagonal_scale", &ksp->viewerDScale, &ksp->formatDScale, &ksp->viewDScale));
387: PetscCall(KSPGetDiagonalScale(ksp, &flag));
388: PetscCall(PetscOptionsBool("-ksp_diagonal_scale", "Diagonal scale matrix before building preconditioner", "KSPSetDiagonalScale", flag, &flag, &flg));
389: if (flg) PetscCall(KSPSetDiagonalScale(ksp, flag));
390: PetscCall(KSPGetDiagonalScaleFix(ksp, &flag));
391: PetscCall(PetscOptionsBool("-ksp_diagonal_scale_fix", "Fix diagonally scaled matrix after solve", "KSPSetDiagonalScaleFix", flag, &flag, &flg));
392: if (flg) PetscCall(KSPSetDiagonalScaleFix(ksp, flag));
393: nmax = ksp->nmax;
394: PetscCall(PetscOptionsDeprecated("-ksp_matsolve_block_size", "-ksp_matsolve_batch_size", "3.15", NULL));
395: PetscCall(PetscOptionsInt("-ksp_matsolve_batch_size", "Maximum number of columns treated simultaneously", "KSPSetMatSolveBatchSize", nmax, &nmax, &flg));
396: if (flg) PetscCall(KSPSetMatSolveBatchSize(ksp, nmax));
397: goto skipoptions;
398: }
400: PetscCall(PetscOptionsInt("-ksp_max_it", "Maximum number of iterations", "KSPSetTolerances", ksp->max_it, &ksp->max_it, NULL));
401: PetscCall(PetscOptionsReal("-ksp_rtol", "Relative decrease in residual norm", "KSPSetTolerances", ksp->rtol, &ksp->rtol, NULL));
402: PetscCall(PetscOptionsReal("-ksp_atol", "Absolute value of residual norm", "KSPSetTolerances", ksp->abstol, &ksp->abstol, NULL));
403: PetscCall(PetscOptionsReal("-ksp_divtol", "Residual norm increase cause divergence", "KSPSetTolerances", ksp->divtol, &ksp->divtol, NULL));
405: PetscCall(PetscOptionsBool("-ksp_converged_use_initial_residual_norm", "Use initial residual norm for computing relative convergence", "KSPConvergedDefaultSetUIRNorm", PETSC_FALSE, &flag, &set));
406: if (set && flag) PetscCall(KSPConvergedDefaultSetUIRNorm(ksp));
407: PetscCall(PetscOptionsBool("-ksp_converged_use_min_initial_residual_norm", "Use minimum of initial residual norm and b for computing relative convergence", "KSPConvergedDefaultSetUMIRNorm", PETSC_FALSE, &flag, &set));
408: if (set && flag) PetscCall(KSPConvergedDefaultSetUMIRNorm(ksp));
409: PetscCall(PetscOptionsBool("-ksp_converged_maxits", "Declare convergence if the maximum number of iterations is reached", "KSPConvergedDefaultSetConvergedMaxits", PETSC_FALSE, &flag, &set));
410: if (set) PetscCall(KSPConvergedDefaultSetConvergedMaxits(ksp, flag));
411: PetscCall(KSPGetConvergedNegativeCurvature(ksp, &flag));
412: PetscCall(PetscOptionsBool("-ksp_converged_neg_curve", "Declare convergence if negative curvature is detected", "KSPConvergedNegativeCurvature", flag, &flag, &set));
413: if (set) PetscCall(KSPSetConvergedNegativeCurvature(ksp, flag));
414: PetscCall(KSPGetReusePreconditioner(ksp, &reuse));
415: PetscCall(PetscOptionsBool("-ksp_reuse_preconditioner", "Use initial preconditioner and don't ever compute a new one", "KSPReusePreconditioner", reuse, &reuse, NULL));
416: PetscCall(KSPSetReusePreconditioner(ksp, reuse));
418: PetscCall(PetscOptionsBool("-ksp_knoll", "Use preconditioner applied to b for initial guess", "KSPSetInitialGuessKnoll", ksp->guess_knoll, &ksp->guess_knoll, NULL));
419: PetscCall(PetscOptionsBool("-ksp_error_if_not_converged", "Generate error if solver does not converge", "KSPSetErrorIfNotConverged", ksp->errorifnotconverged, &ksp->errorifnotconverged, NULL));
420: PetscCall(PetscOptionsFList("-ksp_guess_type", "Initial guess in Krylov method", NULL, KSPGuessList, NULL, guesstype, 256, &flg));
421: if (flg) {
422: PetscCall(KSPGetGuess(ksp, &ksp->guess));
423: PetscCall(KSPGuessSetType(ksp->guess, guesstype));
424: PetscCall(KSPGuessSetFromOptions(ksp->guess));
425: } else { /* old option for KSP */
426: nmax = 2;
427: PetscCall(PetscOptionsIntArray("-ksp_fischer_guess", "Use Paul Fischer's algorithm or its variants for initial guess", "KSPSetUseFischerGuess", model, &nmax, &flag));
428: if (flag) {
429: PetscCheck(nmax == 2, comm, PETSC_ERR_ARG_OUTOFRANGE, "Must pass in model,size as arguments");
430: PetscCall(KSPSetUseFischerGuess(ksp, model[0], model[1]));
431: }
432: }
434: PetscCall(PetscOptionsEList("-ksp_convergence_test", "Convergence test", "KSPSetConvergenceTest", convtests, 3, "default", &indx, &flg));
435: if (flg) {
436: switch (indx) {
437: case 0:
438: PetscCall(KSPConvergedDefaultCreate(&ctx));
439: PetscCall(KSPSetConvergenceTest(ksp, KSPConvergedDefault, ctx, KSPConvergedDefaultDestroy));
440: break;
441: case 1:
442: PetscCall(KSPSetConvergenceTest(ksp, KSPConvergedSkip, NULL, NULL));
443: break;
444: case 2:
445: PetscCall(KSPConvergedDefaultCreate(&ctx));
446: PetscCall(KSPSetConvergenceTest(ksp, KSPLSQRConvergedDefault, ctx, KSPConvergedDefaultDestroy));
447: break;
448: }
449: }
451: PetscCall(KSPSetUpNorms_Private(ksp, PETSC_FALSE, &normtype, NULL));
452: PetscCall(PetscOptionsEnum("-ksp_norm_type", "KSP Norm type", "KSPSetNormType", KSPNormTypes, (PetscEnum)normtype, (PetscEnum *)&normtype, &flg));
453: if (flg) PetscCall(KSPSetNormType(ksp, normtype));
455: PetscCall(PetscOptionsInt("-ksp_check_norm_iteration", "First iteration to compute residual norm", "KSPSetCheckNormIteration", ksp->chknorm, &ksp->chknorm, NULL));
457: PetscCall(PetscOptionsBool("-ksp_lag_norm", "Lag the calculation of the residual norm", "KSPSetLagNorm", ksp->lagnorm, &flag, &flg));
458: if (flg) PetscCall(KSPSetLagNorm(ksp, flag));
460: PetscCall(KSPGetDiagonalScale(ksp, &flag));
461: PetscCall(PetscOptionsBool("-ksp_diagonal_scale", "Diagonal scale matrix before building preconditioner", "KSPSetDiagonalScale", flag, &flag, &flg));
462: if (flg) PetscCall(KSPSetDiagonalScale(ksp, flag));
463: PetscCall(KSPGetDiagonalScaleFix(ksp, &flag));
464: PetscCall(PetscOptionsBool("-ksp_diagonal_scale_fix", "Fix diagonally scaled matrix after solve", "KSPSetDiagonalScaleFix", flag, &flag, &flg));
465: if (flg) PetscCall(KSPSetDiagonalScaleFix(ksp, flag));
467: PetscCall(PetscOptionsBool("-ksp_constant_null_space", "Add constant null space to Krylov solver matrix", "MatSetNullSpace", PETSC_FALSE, &flg, &set));
468: if (set && flg) {
469: MatNullSpace nsp;
470: Mat Amat = NULL;
472: PetscCall(MatNullSpaceCreate(comm, PETSC_TRUE, 0, NULL, &nsp));
473: if (ksp->pc) PetscCall(PCGetOperators(ksp->pc, &Amat, NULL));
474: PetscCheck(Amat, comm, PETSC_ERR_ARG_WRONGSTATE, "Cannot set nullspace, matrix has not yet been provided");
475: PetscCall(MatSetNullSpace(Amat, nsp));
476: PetscCall(MatNullSpaceDestroy(&nsp));
477: }
479: flg = PETSC_FALSE;
480: if (ksp->pc) {
481: PetscCall(PetscObjectTypeCompare((PetscObject)ksp->pc, PCKSP, &flg));
482: if (!flg) PetscCall(PetscObjectTypeCompare((PetscObject)ksp->pc, PCBJACOBI, &flg));
483: if (!flg) PetscCall(PetscObjectTypeCompare((PetscObject)ksp->pc, PCDEFLATION, &flg));
484: }
486: if (flg) {
487: /* Using dynamic tolerance in preconditioner */
488: PetscCall(PetscOptionsString("-sub_ksp_dynamic_tolerance", "Use dynamic tolerance for inner PC", "KSPMonitorDynamicTolerance", "stdout", monfilename, sizeof(monfilename), &flg));
489: if (flg) {
490: void *scale;
491: PetscReal coeff = 1.0;
493: PetscCall(KSPMonitorDynamicToleranceCreate(&scale));
494: PetscCall(PetscOptionsReal("-sub_ksp_dynamic_tolerance", "Coefficient of dynamic tolerance for inner PC", "KSPMonitorDynamicTolerance", coeff, &coeff, &flg));
495: if (flg) PetscCall(KSPMonitorDynamicToleranceSetCoefficient(scale, coeff));
496: PetscCall(KSPMonitorSet(ksp, KSPMonitorDynamicTolerance, scale, KSPMonitorDynamicToleranceDestroy));
497: }
498: }
500: /*
501: Calls Python function
502: */
503: PetscCall(PetscOptionsString("-ksp_monitor_python", "Use Python function", "KSPMonitorSet", NULL, monfilename, sizeof(monfilename), &flg));
504: if (flg) PetscCall(PetscPythonMonitorSet((PetscObject)ksp, monfilename));
505: /*
506: Graphically plots preconditioned residual norm and range of residual element values
507: */
508: PetscCall(PetscOptionsBool("-ksp_monitor_lg_range", "Monitor graphically range of preconditioned residual norm", "KSPMonitorSet", PETSC_FALSE, &flg, &set));
509: if (set && flg) {
510: PetscViewer ctx;
512: PetscCall(PetscViewerDrawOpen(comm, NULL, NULL, PETSC_DECIDE, PETSC_DECIDE, 400, 300, &ctx));
513: PetscCall(KSPMonitorSet(ksp, KSPMonitorLGRange, ctx, (PetscErrorCode(*)(void **))PetscViewerDestroy));
514: }
515: /* TODO Do these show up in help? */
516: PetscCall(PetscOptionsHasName(((PetscObject)ksp)->options, prefix, "-ksp_converged_rate", &flg));
517: if (flg) {
518: const char *RateTypes[] = {"default", "residual", "error", "PetscRateType", "RATE_", NULL};
519: PetscEnum rtype = (PetscEnum)1;
521: PetscCall(PetscOptionsGetEnum(((PetscObject)ksp)->options, prefix, "-ksp_converged_rate_type", RateTypes, &rtype, &flg));
522: if (rtype == (PetscEnum)0 || rtype == (PetscEnum)1) PetscCall(KSPSetResidualHistory(ksp, NULL, PETSC_DETERMINE, PETSC_TRUE));
523: if (rtype == (PetscEnum)0 || rtype == (PetscEnum)2) PetscCall(KSPSetErrorHistory(ksp, NULL, PETSC_DETERMINE, PETSC_TRUE));
524: }
525: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view", &ksp->viewer, &ksp->format, &ksp->view));
526: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_pre", &ksp->viewerPre, &ksp->formatPre, &ksp->viewPre));
528: flg = PETSC_FALSE;
529: PetscCall(PetscOptionsBool("-ksp_converged_reason_view_cancel", "Cancel all the converged reason view functions set using KSPConvergedReasonViewSet", "KSPConvergedReasonViewCancel", PETSC_FALSE, &flg, &set));
530: if (set && flg) PetscCall(KSPConvergedReasonViewCancel(ksp));
531: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_converged_rate", &ksp->viewerRate, &ksp->formatRate, &ksp->viewRate));
532: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_mat", &ksp->viewerMat, &ksp->formatMat, &ksp->viewMat));
533: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_pmat", &ksp->viewerPMat, &ksp->formatPMat, &ksp->viewPMat));
534: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_rhs", &ksp->viewerRhs, &ksp->formatRhs, &ksp->viewRhs));
535: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_solution", &ksp->viewerSol, &ksp->formatSol, &ksp->viewSol));
536: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_mat_explicit", &ksp->viewerMatExp, &ksp->formatMatExp, &ksp->viewMatExp));
537: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_eigenvalues", &ksp->viewerEV, &ksp->formatEV, &ksp->viewEV));
538: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_singularvalues", &ksp->viewerSV, &ksp->formatSV, &ksp->viewSV));
539: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_eigenvalues_explicit", &ksp->viewerEVExp, &ksp->formatEVExp, &ksp->viewEVExp));
540: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_final_residual", &ksp->viewerFinalRes, &ksp->formatFinalRes, &ksp->viewFinalRes));
541: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_preconditioned_operator_explicit", &ksp->viewerPOpExp, &ksp->formatPOpExp, &ksp->viewPOpExp));
542: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_view_diagonal_scale", &ksp->viewerDScale, &ksp->formatDScale, &ksp->viewDScale));
544: /* Deprecated options */
545: if (!ksp->viewEV) {
546: PetscCall(PetscOptionsDeprecated("-ksp_compute_eigenvalues", NULL, "3.9", "Use -ksp_view_eigenvalues"));
547: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_compute_eigenvalues", &ksp->viewerEV, &ksp->formatEV, &ksp->viewEV));
548: }
549: if (!ksp->viewEV) {
550: PetscCall(PetscOptionsDeprecated("-ksp_plot_eigenvalues", NULL, "3.9", "Use -ksp_view_eigenvalues draw"));
551: PetscCall(PetscOptionsName("-ksp_plot_eigenvalues", "[deprecated since PETSc 3.9; use -ksp_view_eigenvalues draw]", "KSPView", &ksp->viewEV));
552: if (ksp->viewEV) {
553: ksp->formatEV = PETSC_VIEWER_DEFAULT;
554: ksp->viewerEV = PETSC_VIEWER_DRAW_(comm);
555: PetscCall(PetscObjectReference((PetscObject)ksp->viewerEV));
556: }
557: }
558: if (!ksp->viewEV) {
559: PetscCall(PetscOptionsDeprecated("-ksp_plot_eigencontours", NULL, "3.9", "Use -ksp_view_eigenvalues draw::draw_contour"));
560: PetscCall(PetscOptionsName("-ksp_plot_eigencontours", "[deprecated since PETSc 3.9; use -ksp_view_eigenvalues draw::draw_contour]", "KSPView", &ksp->viewEV));
561: if (ksp->viewEV) {
562: ksp->formatEV = PETSC_VIEWER_DRAW_CONTOUR;
563: ksp->viewerEV = PETSC_VIEWER_DRAW_(comm);
564: PetscCall(PetscObjectReference((PetscObject)ksp->viewerEV));
565: }
566: }
567: if (!ksp->viewEVExp) {
568: PetscCall(PetscOptionsDeprecated("-ksp_compute_eigenvalues_explicitly", NULL, "3.9", "Use -ksp_view_eigenvalues_explicit"));
569: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_compute_eigenvalues_explicitly", &ksp->viewerEVExp, &ksp->formatEVExp, &ksp->viewEVExp));
570: }
571: if (!ksp->viewEVExp) {
572: PetscCall(PetscOptionsDeprecated("-ksp_plot_eigenvalues_explicitly", NULL, "3.9", "Use -ksp_view_eigenvalues_explicit draw"));
573: PetscCall(PetscOptionsName("-ksp_plot_eigenvalues_explicitly", "[deprecated since PETSc 3.9; use -ksp_view_eigenvalues_explicit draw]", "KSPView", &ksp->viewEVExp));
574: if (ksp->viewEVExp) {
575: ksp->formatEVExp = PETSC_VIEWER_DEFAULT;
576: ksp->viewerEVExp = PETSC_VIEWER_DRAW_(comm);
577: PetscCall(PetscObjectReference((PetscObject)ksp->viewerEVExp));
578: }
579: }
580: if (!ksp->viewSV) {
581: PetscCall(PetscOptionsDeprecated("-ksp_compute_singularvalues", NULL, "3.9", "Use -ksp_view_singularvalues"));
582: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_compute_singularvalues", &ksp->viewerSV, &ksp->formatSV, &ksp->viewSV));
583: }
584: if (!ksp->viewFinalRes) {
585: PetscCall(PetscOptionsDeprecated("-ksp_final_residual", NULL, "3.9", "Use -ksp_view_final_residual"));
586: PetscCall(PetscOptionsGetViewer(comm, ((PetscObject)ksp)->options, prefix, "-ksp_final_residual", &ksp->viewerFinalRes, &ksp->formatFinalRes, &ksp->viewFinalRes));
587: }
589: #if defined(PETSC_HAVE_SAWS)
590: /*
591: Publish convergence information using AMS
592: */
593: PetscCall(PetscOptionsBool("-ksp_monitor_saws", "Publish KSP progress using SAWs", "KSPMonitorSet", PETSC_FALSE, &flg, &set));
594: if (set && flg) {
595: void *ctx;
596: PetscCall(KSPMonitorSAWsCreate(ksp, &ctx));
597: PetscCall(KSPMonitorSet(ksp, KSPMonitorSAWs, ctx, KSPMonitorSAWsDestroy));
598: PetscCall(KSPSetComputeSingularValues(ksp, PETSC_TRUE));
599: }
600: #endif
602: /* -----------------------------------------------------------------------*/
603: PetscCall(KSPSetUpNorms_Private(ksp, PETSC_FALSE, NULL, &pcside));
604: PetscCall(PetscOptionsEnum("-ksp_pc_side", "KSP preconditioner side", "KSPSetPCSide", PCSides, (PetscEnum)pcside, (PetscEnum *)&pcside, &flg));
605: if (flg) PetscCall(KSPSetPCSide(ksp, pcside));
607: if (ksp->viewSV || ksp->viewEV) PetscCall(KSPSetComputeSingularValues(ksp, PETSC_TRUE));
609: #if defined(PETSC_HAVE_SAWS)
610: {
611: PetscBool set;
612: flg = PETSC_FALSE;
613: PetscCall(PetscOptionsBool("-ksp_saws_block", "Block for SAWs at end of KSPSolve", "PetscObjectSAWsBlock", ((PetscObject)ksp)->amspublishblock, &flg, &set));
614: if (set) PetscCall(PetscObjectSAWsSetBlock((PetscObject)ksp, flg));
615: }
616: #endif
618: nmax = ksp->nmax;
619: PetscCall(PetscOptionsDeprecated("-ksp_matsolve_block_size", "-ksp_matsolve_batch_size", "3.15", NULL));
620: PetscCall(PetscOptionsInt("-ksp_matsolve_batch_size", "Maximum number of columns treated simultaneously", "KSPSetMatSolveBatchSize", nmax, &nmax, &flg));
621: if (flg) PetscCall(KSPSetMatSolveBatchSize(ksp, nmax));
623: flg = PETSC_FALSE;
624: PetscCall(PetscOptionsBool("-ksp_use_explicittranspose", "Explicitly transpose the system in KSPSolveTranspose", "KSPSetUseExplicitTranspose", ksp->transpose.use_explicittranspose, &flg, &set));
625: if (set) PetscCall(KSPSetUseExplicitTranspose(ksp, flg));
627: PetscTryTypeMethod(ksp, setfromoptions, PetscOptionsObject);
628: skipoptions:
629: /* process any options handlers added with PetscObjectAddOptionsHandler() */
630: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)ksp, PetscOptionsObject));
631: PetscOptionsEnd();
632: ksp->setfromoptionscalled++;
633: PetscFunctionReturn(PETSC_SUCCESS);
634: }
636: /*@
637: KSPResetFromOptions - Sets `KSP` parameters from user options ONLY if the `KSP` was previously set from options
639: Collective
641: Input Parameter:
642: . ksp - the `KSP` context
644: Level: advanced
646: .seealso: [](ch_ksp), `KSPSetFromOptions()`, `KSPSetOptionsPrefix()`
647: @*/
648: PetscErrorCode KSPResetFromOptions(KSP ksp)
649: {
650: PetscFunctionBegin;
651: if (ksp->setfromoptionscalled) PetscCall(KSPSetFromOptions(ksp));
652: PetscFunctionReturn(PETSC_SUCCESS);
653: }