1 #include <../src/ksp/ksp/impls/gmres/gmresimpl.h> /*I "petscksp.h" I*/
2
3 /*@C
4 KSPGMRESSetOrthogonalization - Sets the orthogonalization routine used by `KSPGMRES` and `KSPFGMRES`.
5
6 Logically Collective
7
8 Input Parameters:
9 + ksp - iterative context obtained from `KSPCreate()`
10 - fcn - orthogonalization function
11
12 Calling sequence of `fcn`:
13 + ksp - the solver context
14 - it - the current iteration
15
16 Options Database Keys:
17 + -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
18 - -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()
19
20 Level: intermediate
21
22 Notes:
23 Two orthogonalization routines are predefined, `KSPGMRESModifiedGramSchmidtOrthogonalization()` and the default
24 `KSPGMRESClassicalGramSchmidtOrthogonalization()`.
25
26 Use `KSPGMRESSetCGSRefinementType()` to determine if iterative refinement is used to increase stability.
27
28 .seealso: [](ch_ksp), `KSPGMRESSetRestart()`, `KSPGMRESSetPreAllocateVectors()`,
29 `KSPGMRESSetCGSRefinementType()`, `KSPGMRESModifiedGramSchmidtOrthogonalization()`,
30 `KSPGMRESClassicalGramSchmidtOrthogonalization()`, `KSPGMRESGetCGSRefinementType()`
31 @*/
KSPGMRESSetOrthogonalization(KSP ksp,PetscErrorCode (* fcn)(KSP ksp,PetscInt it))32 PetscErrorCode KSPGMRESSetOrthogonalization(KSP ksp, PetscErrorCode (*fcn)(KSP ksp, PetscInt it))
33 {
34 PetscFunctionBegin;
35 PetscValidHeaderSpecific(ksp, KSP_CLASSID, 1);
36 PetscTryMethod(ksp, "KSPGMRESSetOrthogonalization_C", (KSP, PetscErrorCode (*)(KSP, PetscInt)), (ksp, fcn));
37 PetscFunctionReturn(PETSC_SUCCESS);
38 }
39
40 /*@C
41 KSPGMRESGetOrthogonalization - Gets the orthogonalization routine used by `KSPGMRES` and `KSPFGMRES`.
42
43 Not Collective
44
45 Input Parameter:
46 . ksp - iterative context obtained from `KSPCreate()`
47
48 Output Parameter:
49 . fcn - orthogonalization function
50
51 Calling sequence of `fcn`:
52 + ksp - the solver context
53 - it - the current iteration
54
55 Level: intermediate
56
57 Notes:
58 Two orthogonalization routines are predefined, `KSPGMRESModifiedGramSchmidtOrthogonalization()`, and the default
59 `KSPGMRESClassicalGramSchmidtOrthogonalization()`
60
61 Use `KSPGMRESSetCGSRefinementType()` to determine if iterative refinement is used to increase stability.
62
63 .seealso: [](ch_ksp), `KSPGMRESSetRestart()`, `KSPGMRESSetPreAllocateVectors()`, `KSPGMRESSetCGSRefinementType()`, `KSPGMRESSetOrthogonalization()`,
64 `KSPGMRESModifiedGramSchmidtOrthogonalization()`, `KSPGMRESClassicalGramSchmidtOrthogonalization()`, `KSPGMRESGetCGSRefinementType()`
65 @*/
KSPGMRESGetOrthogonalization(KSP ksp,PetscErrorCode (** fcn)(KSP ksp,PetscInt it))66 PetscErrorCode KSPGMRESGetOrthogonalization(KSP ksp, PetscErrorCode (**fcn)(KSP ksp, PetscInt it))
67 {
68 PetscFunctionBegin;
69 PetscValidHeaderSpecific(ksp, KSP_CLASSID, 1);
70 PetscUseMethod(ksp, "KSPGMRESGetOrthogonalization_C", (KSP, PetscErrorCode (**)(KSP, PetscInt)), (ksp, fcn));
71 PetscFunctionReturn(PETSC_SUCCESS);
72 }
73