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