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.
 27:    Options Database Keys:
 29: +  -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
 30: -  -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()
 32:    Level: intermediate
 34: .seealso: KSPGMRESSetRestart(), KSPGMRESSetPreAllocateVectors(), KSPGMRESSetCGSRefinementType(), KSPGMRESSetOrthogonalization(),
 35:           KSPGMRESModifiedGramSchmidtOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESGetCGSRefinementType()
 36: @*/
 37: PetscErrorCode  KSPGMRESSetOrthogonalization(KSP ksp,PetscErrorCode (*fcn)(KSP,PetscInt))
 38: {
 43:   PetscTryMethod(ksp,"KSPGMRESSetOrthogonalization_C",(KSP,PetscErrorCode (*)(KSP,PetscInt)),(ksp,fcn));
 44:   return(0);
 45: }
 47: /*@C
 48:    KSPGMRESGetOrthogonalization - Gets the orthogonalization routine used by GMRES and FGMRES.
 50:    Not Collective
 52:    Input Parameter:
 53: .  ksp - iterative context obtained from KSPCreate
 55:    Output Parameter:
 56: .  fcn - orthogonalization function
 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
 63:    Notes:
 64:    Two orthogonalization routines are predefined, including
 66:    KSPGMRESModifiedGramSchmidtOrthogonalization()
 68:    KSPGMRESClassicalGramSchmidtOrthogonalization() - Default. Use KSPGMRESSetCGSRefinementType() to determine if
 69:      iterative refinement is used to increase stability.
 72:    Options Database Keys:
 74: +  -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
 75: -  -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()
 77:    Level: intermediate
 79: .seealso: KSPGMRESSetRestart(), KSPGMRESSetPreAllocateVectors(), KSPGMRESSetCGSRefinementType(), KSPGMRESSetOrthogonalization(),
 80:           KSPGMRESModifiedGramSchmidtOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESGetCGSRefinementType()
 81: @*/
 82: PetscErrorCode  KSPGMRESGetOrthogonalization(KSP ksp,PetscErrorCode (**fcn)(KSP,PetscInt))
 83: {
 88:   PetscUseMethod(ksp,"KSPGMRESGetOrthogonalization_C",(KSP,PetscErrorCode (**)(KSP,PetscInt)),(ksp,fcn));
 89:   return(0);
 90: }