Actual source code: gmres2.c

petsc-3.3-p7 2013-05-11
  2: #include <../src/ksp/ksp/impls/gmres/gmresimpl.h>       /*I  "petscksp.h"  I*/

  6: /*@C
  7:    KSPGMRESSetOrthogonalization - Sets the orthogonalization routine used by GMRES and FGMRES.

  9:    Logically Collective on KSP

 11:    Input Parameters:
 12: +  ksp - iterative context obtained from KSPCreate
 13: -  fcn - orthogonalization function

 15:    Calling Sequence of function:
 16: $   errorcode = int fcn(KSP ksp,int it);
 17: $   it is one minus the number of GMRES iterations since last restart;
 18: $    i.e. the size of Krylov space minus one

 20:    Notes:
 21:    Two orthogonalization routines are predefined, including

 23:    KSPGMRESModifiedGramSchmidtOrthogonalization()

 25:    KSPGMRESClassicalGramSchmidtOrthogonalization() - Default. Use KSPGMRESSetCGSRefinementType() to determine if 
 26:      iterative refinement is used to increase stability. 


 29:    Options Database Keys:

 31: +  -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
 32: -  -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()

 34:    Level: intermediate

 36: .keywords: KSP, GMRES, set, orthogonalization, Gram-Schmidt, iterative refinement

 38: .seealso: KSPGMRESSetRestart(), KSPGMRESSetPreAllocateVectors(), KSPGMRESSetCGSRefinementType(), KSPGMRESSetOrthogonalization(),
 39:           KSPGMRESModifiedGramSchmidtOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESGetCGSRefinementType()
 40: @*/
 41: PetscErrorCode  KSPGMRESSetOrthogonalization(KSP ksp,PetscErrorCode (*fcn)(KSP,PetscInt))
 42: {

 47:   PetscTryMethod(ksp,"KSPGMRESSetOrthogonalization_C",(KSP,PetscErrorCode (*)(KSP,PetscInt)),(ksp,fcn));
 48:   return(0);
 49: }

 53: /*@C
 54:    KSPGMRESGetOrthogonalization - Gets the orthogonalization routine used by GMRES and FGMRES.

 56:    Not Collective

 58:    Input Parameter:
 59: .  ksp - iterative context obtained from KSPCreate

 61:    Output Parameter:
 62: .  fcn - orthogonalization function

 64:    Calling Sequence of function:
 65: $   errorcode = int fcn(KSP ksp,int it);
 66: $   it is one minus the number of GMRES iterations since last restart;
 67: $    i.e. the size of Krylov space minus one

 69:    Notes:
 70:    Two orthogonalization routines are predefined, including

 72:    KSPGMRESModifiedGramSchmidtOrthogonalization()

 74:    KSPGMRESClassicalGramSchmidtOrthogonalization() - Default. Use KSPGMRESSetCGSRefinementType() to determine if 
 75:      iterative refinement is used to increase stability. 


 78:    Options Database Keys:

 80: +  -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
 81: -  -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()

 83:    Level: intermediate

 85: .keywords: KSP, GMRES, set, orthogonalization, Gram-Schmidt, iterative refinement

 87: .seealso: KSPGMRESSetRestart(), KSPGMRESSetPreAllocateVectors(), KSPGMRESSetCGSRefinementType(), KSPGMRESSetOrthogonalization(),
 88:           KSPGMRESModifiedGramSchmidtOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESGetCGSRefinementType()
 89: @*/
 90: PetscErrorCode  KSPGMRESGetOrthogonalization(KSP ksp,PetscErrorCode (**fcn)(KSP,PetscInt))
 91: {

 96:   PetscUseMethod(ksp,"KSPGMRESGetOrthogonalization_C",(KSP,PetscErrorCode (**)(KSP,PetscInt)),(ksp,fcn));
 97:   return(0);
 98: }