Actual source code: gmres2.c
petsc-3.11.0 2019-03-29
2: #include <../src/ksp/ksp/impls/gmres/gmresimpl.h>
4: /*@C
5: KSPGMRESSetOrthogonalization - Sets the orthogonalization routine used by GMRES and FGMRES.
7: Logically Collective on KSP
9: Input Parameters:
10: + ksp - iterative context obtained from KSPCreate
11: - fcn - orthogonalization function
13: Calling Sequence of function:
14: $ errorcode = int fcn(KSP ksp,int it);
15: $ it is one minus the number of GMRES iterations since last restart;
16: $ i.e. the size of Krylov space minus one
18: Notes:
19: Two orthogonalization routines are predefined, including
21: KSPGMRESModifiedGramSchmidtOrthogonalization()
23: KSPGMRESClassicalGramSchmidtOrthogonalization() - Default. Use KSPGMRESSetCGSRefinementType() to determine if
24: iterative refinement is used to increase stability.
27: Options Database Keys:
29: + -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
30: - -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()
32: Level: intermediate
34: .keywords: KSP, GMRES, set, orthogonalization, Gram-Schmidt, iterative refinement
36: .seealso: KSPGMRESSetRestart(), KSPGMRESSetPreAllocateVectors(), KSPGMRESSetCGSRefinementType(), KSPGMRESSetOrthogonalization(),
37: KSPGMRESModifiedGramSchmidtOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESGetCGSRefinementType()
38: @*/
39: PetscErrorCode KSPGMRESSetOrthogonalization(KSP ksp,PetscErrorCode (*fcn)(KSP,PetscInt))
40: {
45: PetscTryMethod(ksp,"KSPGMRESSetOrthogonalization_C",(KSP,PetscErrorCode (*)(KSP,PetscInt)),(ksp,fcn));
46: return(0);
47: }
49: /*@C
50: KSPGMRESGetOrthogonalization - Gets the orthogonalization routine used by GMRES and FGMRES.
52: Not Collective
54: Input Parameter:
55: . ksp - iterative context obtained from KSPCreate
57: Output Parameter:
58: . fcn - orthogonalization function
60: Calling Sequence of function:
61: $ errorcode = int fcn(KSP ksp,int it);
62: $ it is one minus the number of GMRES iterations since last restart;
63: $ i.e. the size of Krylov space minus one
65: Notes:
66: Two orthogonalization routines are predefined, including
68: KSPGMRESModifiedGramSchmidtOrthogonalization()
70: KSPGMRESClassicalGramSchmidtOrthogonalization() - Default. Use KSPGMRESSetCGSRefinementType() to determine if
71: iterative refinement is used to increase stability.
74: Options Database Keys:
76: + -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
77: - -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()
79: Level: intermediate
81: .keywords: KSP, GMRES, set, orthogonalization, Gram-Schmidt, iterative refinement
83: .seealso: KSPGMRESSetRestart(), KSPGMRESSetPreAllocateVectors(), KSPGMRESSetCGSRefinementType(), KSPGMRESSetOrthogonalization(),
84: KSPGMRESModifiedGramSchmidtOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESGetCGSRefinementType()
85: @*/
86: PetscErrorCode KSPGMRESGetOrthogonalization(KSP ksp,PetscErrorCode (**fcn)(KSP,PetscInt))
87: {
92: PetscUseMethod(ksp,"KSPGMRESGetOrthogonalization_C",(KSP,PetscErrorCode (**)(KSP,PetscInt)),(ksp,fcn));
93: return(0);
94: }