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