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: {
 35:   PetscFunctionBegin;
 37:   PetscTryMethod(ksp, "KSPGMRESSetOrthogonalization_C", (KSP, PetscErrorCode(*)(KSP, PetscInt)), (ksp, fcn));
 38:   PetscFunctionReturn(PETSC_SUCCESS);
 39: }

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

 44:    Not Collective

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

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

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

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

 62:    Level: intermediate

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

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

 70: .seealso: [](chapter_ksp), `KSPGMRESSetRestart()`, `KSPGMRESSetPreAllocateVectors()`, `KSPGMRESSetCGSRefinementType()`, `KSPGMRESSetOrthogonalization()`,
 71:           `KSPGMRESModifiedGramSchmidtOrthogonalization()`, `KSPGMRESClassicalGramSchmidtOrthogonalization()`, `KSPGMRESGetCGSRefinementType()`
 72: @*/
 73: PetscErrorCode KSPGMRESGetOrthogonalization(KSP ksp, PetscErrorCode (**fcn)(KSP, PetscInt))
 74: {
 75:   PetscFunctionBegin;
 77:   PetscUseMethod(ksp, "KSPGMRESGetOrthogonalization_C", (KSP, PetscErrorCode(**)(KSP, PetscInt)), (ksp, fcn));
 78:   PetscFunctionReturn(PETSC_SUCCESS);
 79: }