Actual source code: gmres2.c


  2: #include <../src/ksp/ksp/impls/gmres/gmresimpl.h>

  4: /*@C
  5:    KSPGMRESSetOrthogonalization - Sets the orthogonalization routine used by `KSPGMRES` and `KSPFGMRES`.

  7:    Logically Collective

  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:    Options Database Keys:
 19: +  -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
 20: -  -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()

 22:    Level: intermediate

 24:    Notes:
 25:    Two orthogonalization routines are predefined, including `KSPGMRESModifiedGramSchmidtOrthogonalization()` and the default
 26:    `KSPGMRESClassicalGramSchmidtOrthogonalization()`.

 28:    Use `KSPGMRESSetCGSRefinementType()` to determine if iterative refinement is used to increase stability.

 30: .seealso: [](chapter_ksp), `KSPGMRESSetRestart()`, `KSPGMRESSetPreAllocateVectors()`, `KSPGMRESSetCGSRefinementType()`, `KSPGMRESSetOrthogonalization()`,
 31:           `KSPGMRESModifiedGramSchmidtOrthogonalization()`, `KSPGMRESClassicalGramSchmidtOrthogonalization()`, `KSPGMRESGetCGSRefinementType()`
 32: @*/
 33: PetscErrorCode KSPGMRESSetOrthogonalization(KSP ksp, PetscErrorCode (*fcn)(KSP, PetscInt))
 34: {
 36:   PetscTryMethod(ksp, "KSPGMRESSetOrthogonalization_C", (KSP, PetscErrorCode(*)(KSP, PetscInt)), (ksp, fcn));
 37:   return 0;
 38: }

 40: /*@C
 41:    KSPGMRESGetOrthogonalization - Gets the orthogonalization routine used by `KSPGMRES` and `KSPFGMRES`.

 43:    Not Collective

 45:    Input Parameter:
 46: .  ksp - iterative context obtained from `KSPCreate()`

 48:    Output Parameter:
 49: .  fcn - orthogonalization function

 51:    Calling Sequence of function:
 52: .vb
 53:    errorcode = PetscErrorCode fcn(KSP ksp,PetscInt it);
 54:    it is one minus the number of GMRES iterations since last restart; i.e. the size of Krylov space minus one
 55: .ve

 57:    Options Database Keys:
 58: +  -ksp_gmres_classicalgramschmidt - Activates KSPGMRESClassicalGramSchmidtOrthogonalization() (default)
 59: -  -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()

 61:    Level: intermediate

 63:    Notes:
 64:    Two orthogonalization routines are predefined, including `KSPGMRESModifiedGramSchmidtOrthogonalization()`, and the default
 65:    `KSPGMRESClassicalGramSchmidtOrthogonalization()`

 67:    Use `KSPGMRESSetCGSRefinementType()` to determine if iterative refinement is used to increase stability.

 69: .seealso: [](chapter_ksp), `KSPGMRESSetRestart()`, `KSPGMRESSetPreAllocateVectors()`, `KSPGMRESSetCGSRefinementType()`, `KSPGMRESSetOrthogonalization()`,
 70:           `KSPGMRESModifiedGramSchmidtOrthogonalization()`, `KSPGMRESClassicalGramSchmidtOrthogonalization()`, `KSPGMRESGetCGSRefinementType()`
 71: @*/
 72: PetscErrorCode KSPGMRESGetOrthogonalization(KSP ksp, PetscErrorCode (**fcn)(KSP, PetscInt))
 73: {
 75:   PetscUseMethod(ksp, "KSPGMRESGetOrthogonalization_C", (KSP, PetscErrorCode(**)(KSP, PetscInt)), (ksp, fcn));
 76:   return 0;
 77: }