Actual source code: gmres2.c
2: #include <../src/ksp/ksp/impls/gmres/gmresimpl.h>
4: /*@C
5: KSPGMRESSetOrthogonalization - Sets the orthogonalization routine used by `KSPGMRES` and `KSPFGMRES`.
7: Logically Collective
9: Input Parameters:
10: + ksp - iterative context obtained from `KSPCreate()`
11: - fcn - orthogonalization function
13: Calling Sequence of `fcn`:
14: $ PetscErrorCode fcn(KSP ksp, PetscInt it);
15: + KSP - the solver context
16: - it - the current iteration
18: Options Database Keys:
19: + -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
20: - -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()
22: Level: intermediate
24: Notes:
25: Two orthogonalization routines are predefined, including `KSPGMRESModifiedGramSchmidtOrthogonalization()` and the default
26: `KSPGMRESClassicalGramSchmidtOrthogonalization()`.
28: Use `KSPGMRESSetCGSRefinementType()` to determine if iterative refinement is used to increase stability.
30: .seealso: [](ch_ksp), `KSPGMRESSetRestart()`, `KSPGMRESSetPreAllocateVectors()`, `KSPGMRESSetCGSRefinementType()`, `KSPGMRESSetOrthogonalization()`,
31: `KSPGMRESModifiedGramSchmidtOrthogonalization()`, `KSPGMRESClassicalGramSchmidtOrthogonalization()`, `KSPGMRESGetCGSRefinementType()`
32: @*/
33: PetscErrorCode KSPGMRESSetOrthogonalization(KSP ksp, PetscErrorCode (*fcn)(KSP, PetscInt))
34: {
35: PetscFunctionBegin;
37: PetscTryMethod(ksp, "KSPGMRESSetOrthogonalization_C", (KSP, PetscErrorCode(*)(KSP, PetscInt)), (ksp, fcn));
38: PetscFunctionReturn(PETSC_SUCCESS);
39: }
41: /*@C
42: KSPGMRESGetOrthogonalization - Gets the orthogonalization routine used by `KSPGMRES` and `KSPFGMRES`.
44: Not Collective
46: Input Parameter:
47: . ksp - iterative context obtained from `KSPCreate()`
49: Output Parameter:
50: . fcn - orthogonalization function
52: Calling Sequence of `fcn`:
53: .vb
54: PetscErrorCode fcn(KSP ksp, PetscInt it);
55: .ve
56: + KSP - the solver context
57: - it - the current iteration
59: Options Database Keys:
60: + -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
61: - -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()
63: Level: intermediate
65: Notes:
66: Two orthogonalization routines are predefined, including `KSPGMRESModifiedGramSchmidtOrthogonalization()`, and the default
67: `KSPGMRESClassicalGramSchmidtOrthogonalization()`
69: Use `KSPGMRESSetCGSRefinementType()` to determine if iterative refinement is used to increase stability.
71: .seealso: [](ch_ksp), `KSPGMRESSetRestart()`, `KSPGMRESSetPreAllocateVectors()`, `KSPGMRESSetCGSRefinementType()`, `KSPGMRESSetOrthogonalization()`,
72: `KSPGMRESModifiedGramSchmidtOrthogonalization()`, `KSPGMRESClassicalGramSchmidtOrthogonalization()`, `KSPGMRESGetCGSRefinementType()`
73: @*/
74: PetscErrorCode KSPGMRESGetOrthogonalization(KSP ksp, PetscErrorCode (**fcn)(KSP, PetscInt))
75: {
76: PetscFunctionBegin;
78: PetscUseMethod(ksp, "KSPGMRESGetOrthogonalization_C", (KSP, PetscErrorCode(**)(KSP, PetscInt)), (ksp, fcn));
79: PetscFunctionReturn(PETSC_SUCCESS);
80: }