Actual source code: gmres2.c
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 = PetscErrorCode fcn(KSP ksp,PetscInt 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.
26: Options Database Keys:
28: + -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
29: - -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()
31: Level: intermediate
33: .seealso: KSPGMRESSetRestart(), KSPGMRESSetPreAllocateVectors(), KSPGMRESSetCGSRefinementType(), KSPGMRESSetOrthogonalization(),
34: KSPGMRESModifiedGramSchmidtOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESGetCGSRefinementType()
35: @*/
36: PetscErrorCode KSPGMRESSetOrthogonalization(KSP ksp,PetscErrorCode (*fcn)(KSP,PetscInt))
37: {
39: PetscTryMethod(ksp,"KSPGMRESSetOrthogonalization_C",(KSP,PetscErrorCode (*)(KSP,PetscInt)),(ksp,fcn));
40: return 0;
41: }
43: /*@C
44: KSPGMRESGetOrthogonalization - Gets the orthogonalization routine used by GMRES and FGMRES.
46: Not Collective
48: Input Parameter:
49: . ksp - iterative context obtained from KSPCreate
51: Output Parameter:
52: . fcn - orthogonalization function
54: Calling Sequence of function:
55: $ errorcode = PetscErrorCode fcn(KSP ksp,PetscInt it);
56: $ it is one minus the number of GMRES iterations since last restart;
57: $ i.e. the size of Krylov space minus one
59: Notes:
60: Two orthogonalization routines are predefined, including
62: KSPGMRESModifiedGramSchmidtOrthogonalization()
64: KSPGMRESClassicalGramSchmidtOrthogonalization() - Default. Use KSPGMRESSetCGSRefinementType() to determine if
65: iterative refinement is used to increase stability.
67: Options Database Keys:
69: + -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
70: - -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()
72: Level: intermediate
74: .seealso: KSPGMRESSetRestart(), KSPGMRESSetPreAllocateVectors(), KSPGMRESSetCGSRefinementType(), KSPGMRESSetOrthogonalization(),
75: KSPGMRESModifiedGramSchmidtOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESGetCGSRefinementType()
76: @*/
77: PetscErrorCode KSPGMRESGetOrthogonalization(KSP ksp,PetscErrorCode (**fcn)(KSP,PetscInt))
78: {
80: PetscUseMethod(ksp,"KSPGMRESGetOrthogonalization_C",(KSP,PetscErrorCode (**)(KSP,PetscInt)),(ksp,fcn));
81: return 0;
82: }