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 - function to modify the `PC`, see `KSPFlexibleModifyPCFn`
11: . ctx - optional context
12: - destroy - optional context destroy routine
14: Options Database Keys:
15: + -ksp_fgmres_modifypcnochange - do not change the `PC`
16: - -ksp_fgmres_modifypcksp - changes the inner KSP solver tolerances
18: Level: intermediate
20: Note:
21: Several `fcn` routines are predefined, including `KSPFGMRESModifyPCNoChange()` and `KSPFGMRESModifyPCKSP()`
23: .seealso: [](ch_ksp), [](sec_flexibleksp), `KSPFGMRES`, `KSPFlexibleModifyPCFn`, `KSPFlexibleSetModifyPC()`, `KSPFGMRESModifyPCNoChange()`, `KSPFGMRESModifyPCKSP()`
24: @*/
25: PetscErrorCode KSPFGMRESSetModifyPC(KSP ksp, KSPFlexibleModifyPCFn *fcn, void *ctx, PetscCtxDestroyFn *destroy)
26: {
27: PetscFunctionBegin;
29: PetscTryMethod(ksp, "KSPFGMRESSetModifyPC_C", (KSP, KSPFlexibleModifyPCFn *, void *, PetscCtxDestroyFn *), (ksp, fcn, ctx, destroy));
30: PetscFunctionReturn(PETSC_SUCCESS);
31: }
33: /*@C
34: KSPFlexibleSetModifyPC - Sets the routine used by flexible `KSP` methods to modify the preconditioner. [](sec_flexibleksp)
36: Logically Collective
38: Input Parameters:
39: + ksp - iterative context obtained from `KSPCreate()`
40: . fcn - function to modify the `PC`, see `KSPFlexibleModifyPCFn`
41: . ctx - optional context
42: - destroy - optional context destroy routine
44: Level: intermediate
46: Note:
47: Several `fcn` routines are predefined, including `KSPFGMRESModifyPCNoChange()` and `KSPFGMRESModifyPCKSP()`
49: .seealso: [](ch_ksp), [](sec_flexibleksp), `KSPFGMRES`, `KSPFGMRESModifyPCNoChange()`, `KSPFGMRESModifyPCKSP()`
50: @*/
51: PetscErrorCode KSPFlexibleSetModifyPC(KSP ksp, KSPFlexibleModifyPCFn *fcn, void *ctx, PetscCtxDestroyFn *destroy)
52: {
53: PetscFunctionBegin;
55: PetscTryMethod(ksp, "KSPFlexibleSetModifyPC_C", (KSP, KSPFlexibleModifyPCFn *, void *, PetscCtxDestroyFn *), (ksp, fcn, ctx, destroy));
56: PetscFunctionReturn(PETSC_SUCCESS);
57: }
59: /*@C
60: KSPFGMRESModifyPCNoChange - this is the default used by `KSPFMGMRES` - it doesn't change the preconditioner. [](sec_flexibleksp)
62: Input Parameters:
63: + ksp - the ksp context being used.
64: . total_its - the total number of `KSPFGMRES` iterations that have occurred.
65: . loc_its - the number of `KSPFGMRES` iterations since last restart.
66: . res_norm - the current residual norm.
67: - ctx - context variable, unused in this routine
69: Level: intermediate
71: .seealso: [](ch_ksp), [](sec_flexibleksp), `KSPFGMRES`, `KSPFlexibleModifyPCFn`, `KSPFGMRESSetModifyPC()`, `KSPFGMRESModifyPCKSP()`
72: @*/
73: PetscErrorCode KSPFGMRESModifyPCNoChange(KSP ksp, PetscInt total_its, PetscInt loc_its, PetscReal res_norm, void *ctx)
74: {
75: PetscFunctionBegin;
76: PetscFunctionReturn(PETSC_SUCCESS);
77: }
79: /*@C
80: KSPFGMRESModifyPCKSP - modifies the attributes of the `KSPFGMRES` preconditioner, see [](sec_flexibleksp).
82: Input Parameters:
83: + ksp - the ksp context being used.
84: . total_its - the total number of `KSPFGMRES` iterations that have occurred.
85: . loc_its - the number of `KSPFGMRES` iterations since last restart.
86: . res_norm - the current residual norm.
87: - ctx - context, not used in this routine
89: Level: intermediate
91: Note:
92: You can use this as a template for writing a custom modification callback
94: .seealso: [](ch_ksp), [](sec_flexibleksp), `KSPFGMRES`, `KSPFlexibleModifyPCFn`, `KSPFGMRESSetModifyPC()`
95: @*/
96: PetscErrorCode KSPFGMRESModifyPCKSP(KSP ksp, PetscInt total_its, PetscInt loc_its, PetscReal res_norm, void *ctx)
97: {
98: PC pc;
99: PetscInt maxits;
100: KSP sub_ksp;
101: PetscReal rtol, abstol, dtol;
102: PetscBool isksp;
104: PetscFunctionBegin;
105: PetscCall(KSPGetPC(ksp, &pc));
107: PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCKSP, &isksp));
108: if (isksp) {
109: PetscCall(PCKSPGetKSP(pc, &sub_ksp));
111: /* note that at this point you could check the type of KSP with KSPGetType() */
113: /* Now we can use functions such as KSPGMRESSetRestart() or
114: KSPGMRESSetOrthogonalization() or KSPSetTolerances() */
116: PetscCall(KSPGetTolerances(sub_ksp, &rtol, &abstol, &dtol, &maxits));
117: if (!loc_its) rtol = .1;
118: else rtol *= .9;
119: PetscCall(KSPSetTolerances(sub_ksp, rtol, abstol, dtol, maxits));
120: }
121: PetscFunctionReturn(PETSC_SUCCESS);
122: }