Actual source code: modpcf.c
1: #define PETSCKSP_DLL
3: #include petscksp.h
6: /*@C
7: KSPFGMRESSetModifyPC - Sets the routine used by FGMRES to modify the preconditioner.
9: Collective on KSP
11: Input Parameters:
12: + ksp - iterative context obtained from KSPCreate
13: . fcn - modifypc function
14: . ctx - optional contex
15: - d - optional context destroy routine
17: Calling Sequence of function:
18: int fcn(KSP ksp,int total_its,int loc_its,PetscReal res_norm,void*ctx);
20: ksp - the ksp context being used.
21: total_its - the total number of FGMRES iterations that have occurred.
22: loc_its - the number of FGMRES iterations since last restart.
23: res_norm - the current residual norm.
24: ctx - optional context variable
26: Options Database Keys:
27: -ksp_fgmres_modifypcnochange
28: -ksp_fgmres_modifypcksp
30: Level: intermediate
32: Contributed by Allison Baker
34: Notes:
35: Several modifypc routines are predefined, including
36: KSPFGMRESModifyPCNoChange()
37: KSPFGMRESModifyPCKSP()
39: .seealso: KSPFGMRESModifyPCNoChange(), KSPFGMRESModifyPCKSP()
41: @*/
42: PetscErrorCode KSPFGMRESSetModifyPC(KSP ksp,PetscErrorCode (*fcn)(KSP,PetscInt,PetscInt,PetscReal,void*),void* ctx,PetscErrorCode (*d)(void*))
43: {
44: PetscErrorCode ierr,(*f)(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscInt,PetscReal,void*),void*,PetscErrorCode (*)(void*));
48: PetscObjectQueryFunction((PetscObject)ksp,"KSPFGMRESSetModifyPC_C",(void (**)(void))&f);
49: if (f) {
50: (*f)(ksp,fcn,ctx,d);
51: }
52: return(0);
53: }
56: /* The following are different routines used to modify the preconditioner */
60: /*@
62: KSPFGMRESModifyPCNoChange - this is the default used by fgmres - it doesn't change the preconditioner.
64: Input Parameters:
65: + ksp - the ksp context being used.
66: . total_its - the total number of FGMRES iterations that have occurred.
67: . loc_its - the number of FGMRES iterations since last restart.
68: a restart (so number of Krylov directions to be computed)
69: . res_norm - the current residual norm.
70: - dummy - context variable, unused in this routine
72: Level: intermediate
74: Contributed by Allison Baker
76: You can use this as a template!
78: .seealso: KSPFGMRESSetModifyPC(), KSPFGMRESModifyPCKSP()
80: @*/
81: PetscErrorCode KSPFGMRESModifyPCNoChange(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void* dummy)
82: {
85: return(0);
86: }
90: /*@
92: KSPFGMRESModifyPCKSP - modifies the attributes of the
93: GMRES preconditioner. It serves as an example (not as something
94: useful!)
96: Input Parameters:
97: + ksp - the ksp context being used.
98: . total_its - the total number of FGMRES iterations that have occurred.
99: . loc_its - the number of FGMRES iterations since last restart.
100: . res_norm - the current residual norm.
101: - dummy - context, not used here
103: Level: intermediate
105: Contributed by Allison Baker
107: This could be used as a template!
109: .seealso: KSPFGMRESSetModifyPC(), KSPFGMRESModifyPCKSP()
111: @*/
112: PetscErrorCode KSPFGMRESModifyPCKSP(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void *dummy)
113: {
114: PC pc;
116: PetscInt maxits;
117: KSP sub_ksp;
118: PetscReal rtol,abstol,dtol;
119: PetscTruth isksp;
122: KSPGetPC(ksp,&pc);
124: PetscTypeCompare((PetscObject)pc,PCKSP,&isksp);
125: if (isksp) {
126: PCKSPGetKSP(pc,&sub_ksp);
127:
128: /* note that at this point you could check the type of KSP with KSPGetType() */
130: /* Now we can use functions such as KSPGMRESSetRestart() or
131: KSPGMRESSetOrthogonalization() or KSPSetTolerances() */
133: KSPGetTolerances(sub_ksp,&rtol,&abstol,&dtol,&maxits);
134: if (!loc_its) {
135: rtol = .1;
136: } else {
137: rtol *= .9;
138: }
139: KSPSetTolerances(sub_ksp,rtol,abstol,dtol,maxits);
140: }
141: return(0);
142: }