1
2 #include <petsc/private/kspimpl.h> /*I "petscksp.h" I*/
3
4 /*@C
5 KSPFGMRESSetModifyPC - Sets the routine used by FGMRES to modify the preconditioner.
6
7 Logically Collective on ksp
8
9 Input Parameters:
10 + ksp - iterative context obtained from KSPCreate
11 . fcn - modifypc function
12 . ctx - optional contex
13 - d - optional context destroy routine
14
15 Calling Sequence of function:
16 PetscErrorCode fcn(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void*ctx);
17
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
23
24 Options Database Keys:
25 -ksp_fgmres_modifypcnochange
26 -ksp_fgmres_modifypcksp
27
28 Level: intermediate
29
30 Contributed by Allison Baker
31
32 Notes:
33 Several modifypc routines are predefined, including
34 KSPFGMRESModifyPCNoChange()
35 KSPFGMRESModifyPCKSP()
36
37 .seealso: KSPFGMRESModifyPCNoChange(), KSPFGMRESModifyPCKSP()
38
39 @*/
KSPFGMRESSetModifyPC(KSP ksp,PetscErrorCode (* fcn)(KSP,PetscInt,PetscInt,PetscReal,void *),void * ctx,PetscErrorCode (* d)(void *))40 PetscErrorCode KSPFGMRESSetModifyPC(KSP ksp,PetscErrorCode (*fcn)(KSP,PetscInt,PetscInt,PetscReal,void*),void *ctx,PetscErrorCode (*d)(void*))
41 {
42 PetscErrorCode ierr;
43
44 PetscFunctionBegin;
45 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1);
46 ierr = PetscTryMethod(ksp,"KSPFGMRESSetModifyPC_C",(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscInt,PetscReal,void*),void*,PetscErrorCode (*)(void*)),(ksp,fcn,ctx,d));CHKERRQ(ierr);
47 PetscFunctionReturn(0);
48 }
49
50
51 /* The following are different routines used to modify the preconditioner */
52
53 /*@
54
55 KSPFGMRESModifyPCNoChange - this is the default used by fgmres - it doesn't change the preconditioner.
56
57 Input Parameters:
58 + ksp - the ksp context being used.
59 . total_its - the total number of FGMRES iterations that have occurred.
60 . loc_its - the number of FGMRES iterations since last restart.
61 a restart (so number of Krylov directions to be computed)
62 . res_norm - the current residual norm.
63 - dummy - context variable, unused in this routine
64
65 Level: intermediate
66
67 Contributed by Allison Baker
68
69 You can use this as a template!
70
71 .seealso: KSPFGMRESSetModifyPC(), KSPFGMRESModifyPCKSP()
72
73 @*/
KSPFGMRESModifyPCNoChange(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void * dummy)74 PetscErrorCode KSPFGMRESModifyPCNoChange(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void *dummy)
75 {
76 PetscFunctionBegin;
77 PetscFunctionReturn(0);
78 }
79
80 /*@
81
82 KSPFGMRESModifyPCKSP - modifies the attributes of the
83 GMRES preconditioner. It serves as an example (not as something
84 useful!)
85
86 Input Parameters:
87 + ksp - the ksp context being used.
88 . total_its - the total number of FGMRES iterations that have occurred.
89 . loc_its - the number of FGMRES iterations since last restart.
90 . res_norm - the current residual norm.
91 - dummy - context, not used here
92
93 Level: intermediate
94
95 Contributed by Allison Baker
96
97 This could be used as a template!
98
99 .seealso: KSPFGMRESSetModifyPC(), KSPFGMRESModifyPCKSP()
100
101 @*/
KSPFGMRESModifyPCKSP(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void * dummy)102 PetscErrorCode KSPFGMRESModifyPCKSP(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void *dummy)
103 {
104 PC pc;
105 PetscErrorCode ierr;
106 PetscInt maxits;
107 KSP sub_ksp;
108 PetscReal rtol,abstol,dtol;
109 PetscBool isksp;
110
111 PetscFunctionBegin;
112 ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
113
114 ierr = PetscObjectTypeCompare((PetscObject)pc,PCKSP,&isksp);CHKERRQ(ierr);
115 if (isksp) {
116 ierr = PCKSPGetKSP(pc,&sub_ksp);CHKERRQ(ierr);
117
118 /* note that at this point you could check the type of KSP with KSPGetType() */
119
120 /* Now we can use functions such as KSPGMRESSetRestart() or
121 KSPGMRESSetOrthogonalization() or KSPSetTolerances() */
122
123 ierr = KSPGetTolerances(sub_ksp,&rtol,&abstol,&dtol,&maxits);CHKERRQ(ierr);
124 if (!loc_its) rtol = .1;
125 else rtol *= .9;
126 ierr = KSPSetTolerances(sub_ksp,rtol,abstol,dtol,maxits);CHKERRQ(ierr);
127 }
128 PetscFunctionReturn(0);
129 }
130
131
132
133
134