Actual source code: modpcf.c
2: #include <petsc/private/kspimpl.h>
4: /*@C
5: KSPFGMRESSetModifyPC - Sets the routine used by FGMRES to modify the preconditioner.
7: Logically Collective on ksp
9: Input Parameters:
10: + ksp - iterative context obtained from KSPCreate
11: . fcn - modifypc function
12: . ctx - optional context
13: - d - optional context destroy routine
15: Calling Sequence of function:
16: PetscErrorCode fcn(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void*ctx);
18: ksp - the ksp context being used.
19: total_its - the total number of FGMRES iterations that have occurred.
20: loc_its - the number of FGMRES iterations since last restart.
21: res_norm - the current residual norm.
22: ctx - optional context variable
24: Options Database Keys:
25: -ksp_fgmres_modifypcnochange
26: -ksp_fgmres_modifypcksp
28: Level: intermediate
30: Contributed by Allison Baker
32: Notes:
33: Several modifypc routines are predefined, including
34: KSPFGMRESModifyPCNoChange()
35: KSPFGMRESModifyPCKSP()
37: .seealso: KSPFGMRESModifyPCNoChange(), KSPFGMRESModifyPCKSP()
39: @*/
40: PetscErrorCode KSPFGMRESSetModifyPC(KSP ksp,PetscErrorCode (*fcn)(KSP,PetscInt,PetscInt,PetscReal,void*),void *ctx,PetscErrorCode (*d)(void*))
41: {
43: PetscTryMethod(ksp,"KSPFGMRESSetModifyPC_C",(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscInt,PetscReal,void*),void*,PetscErrorCode (*)(void*)),(ksp,fcn,ctx,d));
44: return 0;
45: }
47: /* The following are different routines used to modify the preconditioner */
49: /*@
51: KSPFGMRESModifyPCNoChange - this is the default used by fgmres - it doesn't change the preconditioner.
53: Input Parameters:
54: + ksp - the ksp context being used.
55: . total_its - the total number of FGMRES iterations that have occurred.
56: . loc_its - the number of FGMRES iterations since last restart.
57: a restart (so number of Krylov directions to be computed)
58: . res_norm - the current residual norm.
59: - dummy - context variable, unused in this routine
61: Level: intermediate
63: Contributed by Allison Baker
65: You can use this as a template!
67: .seealso: KSPFGMRESSetModifyPC(), KSPFGMRESModifyPCKSP()
69: @*/
70: PetscErrorCode KSPFGMRESModifyPCNoChange(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void *dummy)
71: {
72: return 0;
73: }
75: /*@
77: KSPFGMRESModifyPCKSP - modifies the attributes of the
78: GMRES preconditioner. It serves as an example (not as something
79: useful!)
81: Input Parameters:
82: + ksp - the ksp context being used.
83: . total_its - the total number of FGMRES iterations that have occurred.
84: . loc_its - the number of FGMRES iterations since last restart.
85: . res_norm - the current residual norm.
86: - dummy - context, not used here
88: Level: intermediate
90: Contributed by Allison Baker
92: This could be used as a template!
94: .seealso: KSPFGMRESSetModifyPC(), KSPFGMRESModifyPCKSP()
96: @*/
97: PetscErrorCode KSPFGMRESModifyPCKSP(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void *dummy)
98: {
99: PC pc;
100: PetscInt maxits;
101: KSP sub_ksp;
102: PetscReal rtol,abstol,dtol;
103: PetscBool isksp;
105: KSPGetPC(ksp,&pc);
107: PetscObjectTypeCompare((PetscObject)pc,PCKSP,&isksp);
108: if (isksp) {
109: 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: KSPGetTolerances(sub_ksp,&rtol,&abstol,&dtol,&maxits);
117: if (!loc_its) rtol = .1;
118: else rtol *= .9;
119: KSPSetTolerances(sub_ksp,rtol,abstol,dtol,maxits);
120: }
121: return 0;
122: }