1
2 #include <../src/ksp/ksp/impls/gmres/gmresimpl.h> /*I "petscksp.h" I*/
3
4 /*@C
5 KSPGMRESSetOrthogonalization - Sets the orthogonalization routine used by GMRES and FGMRES.
6
7 Logically Collective on ksp
8
9 Input Parameters:
10 + ksp - iterative context obtained from KSPCreate
11 - fcn - orthogonalization function
12
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
17
18 Notes:
19 Two orthogonalization routines are predefined, including
20
21 KSPGMRESModifiedGramSchmidtOrthogonalization()
22
23 KSPGMRESClassicalGramSchmidtOrthogonalization() - Default. Use KSPGMRESSetCGSRefinementType() to determine if
24 iterative refinement is used to increase stability.
25
26
27 Options Database Keys:
28
29 + -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
30 - -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()
31
32 Level: intermediate
33
34 .seealso: KSPGMRESSetRestart(), KSPGMRESSetPreAllocateVectors(), KSPGMRESSetCGSRefinementType(), KSPGMRESSetOrthogonalization(),
35 KSPGMRESModifiedGramSchmidtOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESGetCGSRefinementType()
36 @*/
KSPGMRESSetOrthogonalization(KSP ksp,PetscErrorCode (* fcn)(KSP,PetscInt))37 PetscErrorCode KSPGMRESSetOrthogonalization(KSP ksp,PetscErrorCode (*fcn)(KSP,PetscInt))
38 {
39 PetscErrorCode ierr;
40
41 PetscFunctionBegin;
42 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1);
43 ierr = PetscTryMethod(ksp,"KSPGMRESSetOrthogonalization_C",(KSP,PetscErrorCode (*)(KSP,PetscInt)),(ksp,fcn));CHKERRQ(ierr);
44 PetscFunctionReturn(0);
45 }
46
47 /*@C
48 KSPGMRESGetOrthogonalization - Gets the orthogonalization routine used by GMRES and FGMRES.
49
50 Not Collective
51
52 Input Parameter:
53 . ksp - iterative context obtained from KSPCreate
54
55 Output Parameter:
56 . fcn - orthogonalization function
57
58 Calling Sequence of function:
59 $ errorcode = PetscErrorCode fcn(KSP ksp,PetscInt it);
60 $ it is one minus the number of GMRES iterations since last restart;
61 $ i.e. the size of Krylov space minus one
62
63 Notes:
64 Two orthogonalization routines are predefined, including
65
66 KSPGMRESModifiedGramSchmidtOrthogonalization()
67
68 KSPGMRESClassicalGramSchmidtOrthogonalization() - Default. Use KSPGMRESSetCGSRefinementType() to determine if
69 iterative refinement is used to increase stability.
70
71
72 Options Database Keys:
73
74 + -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
75 - -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()
76
77 Level: intermediate
78
79 .seealso: KSPGMRESSetRestart(), KSPGMRESSetPreAllocateVectors(), KSPGMRESSetCGSRefinementType(), KSPGMRESSetOrthogonalization(),
80 KSPGMRESModifiedGramSchmidtOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESGetCGSRefinementType()
81 @*/
KSPGMRESGetOrthogonalization(KSP ksp,PetscErrorCode (** fcn)(KSP,PetscInt))82 PetscErrorCode KSPGMRESGetOrthogonalization(KSP ksp,PetscErrorCode (**fcn)(KSP,PetscInt))
83 {
84 PetscErrorCode ierr;
85
86 PetscFunctionBegin;
87 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1);
88 ierr = PetscUseMethod(ksp,"KSPGMRESGetOrthogonalization_C",(KSP,PetscErrorCode (**)(KSP,PetscInt)),(ksp,fcn));CHKERRQ(ierr);
89 PetscFunctionReturn(0);
90 }
91