Actual source code: pcksp.c
1: #define PETSCKSP_DLL
3: #include private/pcimpl.h
4: #include petscksp.h
6: typedef struct {
7: PetscTruth use_true_matrix; /* use mat rather than pmat in inner linear solve */
8: KSP ksp;
9: PetscInt its; /* total number of iterations KSP uses */
10: } PC_KSP;
14: static PetscErrorCode PCApply_KSP(PC pc,Vec x,Vec y)
15: {
17: PetscInt its;
18: PC_KSP *jac = (PC_KSP*)pc->data;
21: KSPSetInitialGuessNonzero(jac->ksp,pc->nonzero_guess);
22: KSPSolve(jac->ksp,x,y);
23: KSPGetIterationNumber(jac->ksp,&its);
24: jac->its += its;
25: return(0);
26: }
30: static PetscErrorCode PCApplyTranspose_KSP(PC pc,Vec x,Vec y)
31: {
33: PetscInt its;
34: PC_KSP *jac = (PC_KSP*)pc->data;
37: KSPSolveTranspose(jac->ksp,x,y);
38: KSPGetIterationNumber(jac->ksp,&its);
39: jac->its += its;
40: return(0);
41: }
45: static PetscErrorCode PCSetUp_KSP(PC pc)
46: {
48: PC_KSP *jac = (PC_KSP*)pc->data;
49: Mat mat;
50: PetscTruth A;
53: KSPSetFromOptions(jac->ksp);
54: if (jac->use_true_matrix) mat = pc->mat;
55: else mat = pc->pmat;
57: KSPGetOperatorsSet(jac->ksp,&A,PETSC_NULL);
58: if (!A) {
59: KSPSetOperators(jac->ksp,mat,pc->pmat,pc->flag);
60: }
61: KSPSetUp(jac->ksp);
62: return(0);
63: }
65: /* Default destroy, if it has never been setup */
68: static PetscErrorCode PCDestroy_KSP(PC pc)
69: {
70: PC_KSP *jac = (PC_KSP*)pc->data;
74: KSPDestroy(jac->ksp);
75: PetscFree(jac);
76: return(0);
77: }
81: static PetscErrorCode PCView_KSP(PC pc,PetscViewer viewer)
82: {
83: PC_KSP *jac = (PC_KSP*)pc->data;
85: PetscTruth iascii;
88: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
89: if (iascii) {
90: if (jac->use_true_matrix) {
91: PetscViewerASCIIPrintf(viewer,"Using true matrix (not preconditioner matrix) on inner solve\n");
92: }
93: PetscViewerASCIIPrintf(viewer,"KSP and PC on KSP preconditioner follow\n");
94: PetscViewerASCIIPrintf(viewer,"---------------------------------\n");
95: } else {
96: SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for this object",((PetscObject)viewer)->type_name);
97: }
98: PetscViewerASCIIPushTab(viewer);
99: KSPView(jac->ksp,viewer);
100: PetscViewerASCIIPopTab(viewer);
101: if (iascii) {
102: PetscViewerASCIIPrintf(viewer,"---------------------------------\n");
103: }
104: return(0);
105: }
109: static PetscErrorCode PCSetFromOptions_KSP(PC pc)
110: {
112: PetscTruth flg;
115: PetscOptionsHead("KSP preconditioner options");
116: PetscOptionsName("-pc_ksp_true","Use true matrix to define inner linear system, not preconditioner matrix","PCKSPSetUseTrue",&flg);
117: if (flg) {
118: PCKSPSetUseTrue(pc);
119: }
120: PetscOptionsTail();
121: return(0);
122: }
124: /* ----------------------------------------------------------------------------------*/
129: PetscErrorCode PCKSPSetUseTrue_KSP(PC pc)
130: {
131: PC_KSP *jac;
134: jac = (PC_KSP*)pc->data;
135: jac->use_true_matrix = PETSC_TRUE;
136: return(0);
137: }
143: PetscErrorCode PCKSPGetKSP_KSP(PC pc,KSP *ksp)
144: {
145: PC_KSP *jac;
148: jac = (PC_KSP*)pc->data;
149: *ksp = jac->ksp;
150: return(0);
151: }
156: /*@
157: PCKSPSetUseTrue - Sets a flag to indicate that the true matrix (rather than
158: the matrix used to define the preconditioner) is used to compute the
159: residual inside the inner solve.
161: Collective on PC
163: Input Parameters:
164: . pc - the preconditioner context
166: Options Database Key:
167: . -pc_ksp_true - Activates PCKSPSetUseTrue()
169: Note:
170: For the common case in which the preconditioning and linear
171: system matrices are identical, this routine is unnecessary.
173: Level: advanced
175: .keywords: PC, KSP, set, true, local, flag
177: .seealso: PCSetOperators(), PCBJacobiSetUseTrueLocal()
178: @*/
179: PetscErrorCode PCKSPSetUseTrue(PC pc)
180: {
181: PetscErrorCode ierr,(*f)(PC);
185: PetscObjectQueryFunction((PetscObject)pc,"PCKSPSetUseTrue_C",(void (**)(void))&f);
186: if (f) {
187: (*f)(pc);
188: }
189: return(0);
190: }
194: /*@
195: PCKSPGetKSP - Gets the KSP context for a KSP PC.
197: Not Collective but KSP returned is parallel if PC was parallel
199: Input Parameter:
200: . pc - the preconditioner context
202: Output Parameters:
203: . ksp - the PC solver
205: Notes:
206: You must call KSPSetUp() before calling PCKSPGetKSP().
208: Level: advanced
210: .keywords: PC, KSP, get, context
211: @*/
212: PetscErrorCode PCKSPGetKSP(PC pc,KSP *ksp)
213: {
214: PetscErrorCode ierr,(*f)(PC,KSP*);
219: PetscObjectQueryFunction((PetscObject)pc,"PCKSPGetKSP_C",(void (**)(void))&f);
220: if (f) {
221: (*f)(pc,ksp);
222: }
223: return(0);
224: }
226: /* ----------------------------------------------------------------------------------*/
228: /*MC
229: PCKSP - Defines a preconditioner that can consist of any KSP solver.
230: This allows, for example, embedding a Krylov method inside a preconditioner.
232: Options Database Key:
233: . -pc_ksp_true - use the matrix that defines the linear system as the matrix for the
234: inner solver, otherwise by default it uses the matrix used to construct
235: the preconditioner (see PCSetOperators())
237: Level: intermediate
239: Concepts: inner iteration
241: Notes: Using a Krylov method inside another Krylov method can be dangerous (you get divergence or
242: the incorrect answer) unless you use KSPFGMRES as the other Krylov method
245: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
246: PCSHELL, PCCOMPOSITE, PCKSPUseTrue(), PCKSPGetKSP()
248: M*/
253: PetscErrorCode PCCreate_KSP(PC pc)
254: {
256: const char *prefix;
257: PC_KSP *jac;
260: PetscNewLog(pc,PC_KSP,&jac);
261: pc->ops->apply = PCApply_KSP;
262: pc->ops->applytranspose = PCApplyTranspose_KSP;
263: pc->ops->setup = PCSetUp_KSP;
264: pc->ops->destroy = PCDestroy_KSP;
265: pc->ops->setfromoptions = PCSetFromOptions_KSP;
266: pc->ops->view = PCView_KSP;
267: pc->ops->applyrichardson = 0;
269: pc->data = (void*)jac;
270: KSPCreate(((PetscObject)pc)->comm,&jac->ksp);
272: PCGetOptionsPrefix(pc,&prefix);
273: KSPSetOptionsPrefix(jac->ksp,prefix);
274: KSPAppendOptionsPrefix(jac->ksp,"ksp_");
275: jac->use_true_matrix = PETSC_FALSE;
276: jac->its = 0;
278: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCKSPSetUseTrue_C","PCKSPSetUseTrue_KSP",
279: PCKSPSetUseTrue_KSP);
280: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCKSPGetKSP_C","PCKSPGetKSP_KSP",
281: PCKSPGetKSP_KSP);
283: return(0);
284: }