Actual source code: modpcf.c
  1: #include <petsc/private/kspimpl.h>
  3: /*@C
  4:   KSPFGMRESSetModifyPC - Sets the routine used by `KSPFGMRES` to modify the preconditioner. [](sec_flexibleksp)
  6:   Logically Collective
  8:   Input Parameters:
  9: + ksp     - iterative context obtained from `KSPCreate()`
 10: . fcn     - modifypc function
 11: . ctx     - optional context
 12: - destroy - optional context destroy routine
 14:   Calling sequence of `fcn`:
 15: + ksp       - the ksp context being used.
 16: . total_its - the total number of FGMRES iterations that have occurred.
 17: . local_its - the number of FGMRES iterations since last restart.
 18: . res_norm  - the current residual norm.
 19: - ctx       - optional context variable
 21:   Calling sequence of `destroy`:
 22: . ctx - the function context
 24:   Options Database Keys:
 25: + -ksp_fgmres_modifypcnochange - do not change the `PC`
 26: - -ksp_fgmres_modifypcksp      - changes the inner KSP solver tolerances
 28:   Level: intermediate
 30:   Note:
 31:   Several `fcn` routines are predefined, including `KSPFGMRESModifyPCNoChange()`, and  `KSPFGMRESModifyPCKSP()`
 33: .seealso: [](ch_ksp), [](sec_flexibleksp), `KSPFGMRES`, `KSPFGMRESModifyPCNoChange()`, `KSPFGMRESModifyPCKSP()`
 34: @*/
 35: PetscErrorCode KSPFGMRESSetModifyPC(KSP ksp, PetscErrorCode (*fcn)(KSP ksp, PetscInt total_its, PetscInt local_its, PetscReal res_norm, void *ctx), void *ctx, PetscErrorCode (*destroy)(void *ctx))
 36: {
 37:   PetscFunctionBegin;
 39:   PetscTryMethod(ksp, "KSPFGMRESSetModifyPC_C", (KSP, PetscErrorCode (*)(KSP, PetscInt, PetscInt, PetscReal, void *), void *, PetscErrorCode (*)(void *)), (ksp, fcn, ctx, destroy));
 40:   PetscFunctionReturn(PETSC_SUCCESS);
 41: }
 43: /*@
 44:   KSPFGMRESModifyPCNoChange - this is the default used by `KSPFMGMRES` - it doesn't change the preconditioner. [](sec_flexibleksp)
 46:   Input Parameters:
 47: + ksp       - the ksp context being used.
 48: . total_its - the total number of `KSPFGMRES` iterations that have occurred.
 49: . loc_its   - the number of `KSPFGMRES` iterations since last restart.
 50: . res_norm  - the current residual norm.
 51: - ctx       - context variable, unused in this routine
 53:   Level: intermediate
 55: .seealso: [](ch_ksp), [](sec_flexibleksp), `KSPFGMRES`, `KSPFGMRESSetModifyPC()`, `KSPFGMRESModifyPCKSP()`
 56: @*/
 57: PetscErrorCode KSPFGMRESModifyPCNoChange(KSP ksp, PetscInt total_its, PetscInt loc_its, PetscReal res_norm, void *ctx)
 58: {
 59:   PetscFunctionBegin;
 60:   PetscFunctionReturn(PETSC_SUCCESS);
 61: }
 63: /*@
 64:   KSPFGMRESModifyPCKSP - modifies the attributes of the `KSPFGMRES` preconditioner, see [](sec_flexibleksp).
 66:   Input Parameters:
 67: + ksp       - the ksp context being used.
 68: . total_its - the total number of `KSPFGMRES` iterations that have occurred.
 69: . loc_its   - the number of `KSPFGMRES` iterations since last restart.
 70: . res_norm  - the current residual norm.
 71: - ctx       - context, not used here
 73:   Level: intermediate
 75:   Note:
 76:   You can use this as a template for writing a custom modification callback
 78: .seealso: [](ch_ksp), [](sec_flexibleksp), `KSPFGMRES`, `KSPFGMRESSetModifyPC()`
 79: @*/
 80: PetscErrorCode KSPFGMRESModifyPCKSP(KSP ksp, PetscInt total_its, PetscInt loc_its, PetscReal res_norm, void *ctx)
 81: {
 82:   PC        pc;
 83:   PetscInt  maxits;
 84:   KSP       sub_ksp;
 85:   PetscReal rtol, abstol, dtol;
 86:   PetscBool isksp;
 88:   PetscFunctionBegin;
 89:   PetscCall(KSPGetPC(ksp, &pc));
 91:   PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCKSP, &isksp));
 92:   if (isksp) {
 93:     PetscCall(PCKSPGetKSP(pc, &sub_ksp));
 95:     /* note that at this point you could check the type of KSP with KSPGetType() */
 97:     /* Now we can use functions such as KSPGMRESSetRestart() or
 98:       KSPGMRESSetOrthogonalization() or KSPSetTolerances() */
100:     PetscCall(KSPGetTolerances(sub_ksp, &rtol, &abstol, &dtol, &maxits));
101:     if (!loc_its) rtol = .1;
102:     else rtol *= .9;
103:     PetscCall(KSPSetTolerances(sub_ksp, rtol, abstol, dtol, maxits));
104:   }
105:   PetscFunctionReturn(PETSC_SUCCESS);
106: }