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: KSPSolve(jac->ksp,x,y);
22: KSPGetIterationNumber(jac->ksp,&its);
23: jac->its += its;
24: return(0);
25: }
29: static PetscErrorCode PCApplyTranspose_KSP(PC pc,Vec x,Vec y)
30: {
32: PetscInt its;
33: PC_KSP *jac = (PC_KSP*)pc->data;
36: KSPSolveTranspose(jac->ksp,x,y);
37: KSPGetIterationNumber(jac->ksp,&its);
38: jac->its += its;
39: return(0);
40: }
44: static PetscErrorCode PCSetUp_KSP(PC pc)
45: {
47: PC_KSP *jac = (PC_KSP*)pc->data;
48: Mat mat,A;
51: KSPSetFromOptions(jac->ksp);
52: if (jac->use_true_matrix) mat = pc->mat;
53: else mat = pc->pmat;
55: KSPGetOperators(jac->ksp,&A,PETSC_NULL,PETSC_NULL);
56: if (!A) {
57: KSPSetOperators(jac->ksp,mat,pc->pmat,pc->flag);
58: }
59: KSPSetUp(jac->ksp);
60: return(0);
61: }
63: /* Default destroy, if it has never been setup */
66: static PetscErrorCode PCDestroy_KSP(PC pc)
67: {
68: PC_KSP *jac = (PC_KSP*)pc->data;
72: KSPDestroy(jac->ksp);
73: PetscFree(jac);
74: return(0);
75: }
79: static PetscErrorCode PCView_KSP(PC pc,PetscViewer viewer)
80: {
81: PC_KSP *jac = (PC_KSP*)pc->data;
83: PetscTruth iascii;
86: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
87: if (iascii) {
88: if (jac->use_true_matrix) {
89: PetscViewerASCIIPrintf(viewer,"Using true matrix (not preconditioner matrix) on inner solve\n");
90: }
91: PetscViewerASCIIPrintf(viewer,"KSP and PC on KSP preconditioner follow\n");
92: PetscViewerASCIIPrintf(viewer,"---------------------------------\n");
93: } else {
94: SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for this object",((PetscObject)viewer)->type_name);
95: }
96: PetscViewerASCIIPushTab(viewer);
97: KSPView(jac->ksp,viewer);
98: PetscViewerASCIIPopTab(viewer);
99: if (iascii) {
100: PetscViewerASCIIPrintf(viewer,"---------------------------------\n");
101: }
102: return(0);
103: }
107: static PetscErrorCode PCSetFromOptions_KSP(PC pc)
108: {
110: PetscTruth flg;
113: PetscOptionsHead("KSP preconditioner options");
114: PetscOptionsName("-pc_ksp_true","Use true matrix to define inner linear system, not preconditioner matrix","PCKSPSetUseTrue",&flg);
115: if (flg) {
116: PCKSPSetUseTrue(pc);
117: }
118: PetscOptionsTail();
119: return(0);
120: }
122: /* ----------------------------------------------------------------------------------*/
127: PetscErrorCode PETSCKSP_DLLEXPORT PCKSPSetUseTrue_KSP(PC pc)
128: {
129: PC_KSP *jac;
132: jac = (PC_KSP*)pc->data;
133: jac->use_true_matrix = PETSC_TRUE;
134: return(0);
135: }
141: PetscErrorCode PETSCKSP_DLLEXPORT PCKSPGetKSP_KSP(PC pc,KSP *ksp)
142: {
143: PC_KSP *jac;
146: jac = (PC_KSP*)pc->data;
147: *ksp = jac->ksp;
148: return(0);
149: }
154: /*@
155: PCKSPSetUseTrue - Sets a flag to indicate that the true matrix (rather than
156: the matrix used to define the preconditioner) is used to compute the
157: residual inside the inner solve.
159: Collective on PC
161: Input Parameters:
162: . pc - the preconditioner context
164: Options Database Key:
165: . -pc_ksp_true - Activates PCKSPSetUseTrue()
167: Note:
168: For the common case in which the preconditioning and linear
169: system matrices are identical, this routine is unnecessary.
171: Level: advanced
173: .keywords: PC, KSP, set, true, local, flag
175: .seealso: PCSetOperators(), PCBJacobiSetUseTrueLocal()
176: @*/
177: PetscErrorCode PETSCKSP_DLLEXPORT PCKSPSetUseTrue(PC pc)
178: {
179: PetscErrorCode ierr,(*f)(PC);
183: PetscObjectQueryFunction((PetscObject)pc,"PCKSPSetUseTrue_C",(void (**)(void))&f);
184: if (f) {
185: (*f)(pc);
186: }
187: return(0);
188: }
192: /*@
193: PCKSPGetKSP - Gets the KSP context for a KSP PC.
195: Not Collective but KSP returned is parallel if PC was parallel
197: Input Parameter:
198: . pc - the preconditioner context
200: Output Parameters:
201: . ksp - the PC solver
203: Notes:
204: You must call KSPSetUp() before calling PCKSPGetKSP().
206: Level: advanced
208: .keywords: PC, KSP, get, context
209: @*/
210: PetscErrorCode PETSCKSP_DLLEXPORT PCKSPGetKSP(PC pc,KSP *ksp)
211: {
212: PetscErrorCode ierr,(*f)(PC,KSP*);
217: PetscObjectQueryFunction((PetscObject)pc,"PCKSPGetKSP_C",(void (**)(void))&f);
218: if (f) {
219: (*f)(pc,ksp);
220: }
221: return(0);
222: }
224: /* ----------------------------------------------------------------------------------*/
226: /*MC
227: PCKSP - Defines a preconditioner that can consist of any KSP solver.
228: This allows, for example, embedding a Krylov method inside a preconditioner.
230: Options Database Key:
231: . -pc_ksp_true - use the matrix that defines the linear system as the matrix for the
232: inner solver, otherwise by default it uses the matrix used to construct
233: the preconditioner (see PCSetOperators())
235: Level: intermediate
237: Concepts: inner iteration
239: Notes: Using a Krylov method inside another Krylov method can be dangerous (you get divergence or
240: the incorrect answer) unless you use KSPFGMRES as the other Krylov method
243: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
244: PCSHELL, PCCOMPOSITE, PCKSPUseTrue(), PCKSPGetKSP()
246: M*/
251: PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_KSP(PC pc)
252: {
254: const char *prefix;
255: PC_KSP *jac;
258: PetscNew(PC_KSP,&jac);
259: PetscLogObjectMemory(pc,sizeof(PC_KSP));
260: pc->ops->apply = PCApply_KSP;
261: pc->ops->applytranspose = PCApplyTranspose_KSP;
262: pc->ops->setup = PCSetUp_KSP;
263: pc->ops->destroy = PCDestroy_KSP;
264: pc->ops->setfromoptions = PCSetFromOptions_KSP;
265: pc->ops->view = PCView_KSP;
266: pc->ops->applyrichardson = 0;
268: pc->data = (void*)jac;
269: KSPCreate(pc->comm,&jac->ksp);
271: PCGetOptionsPrefix(pc,&prefix);
272: KSPSetOptionsPrefix(jac->ksp,prefix);
273: KSPAppendOptionsPrefix(jac->ksp,"ksp_");
274: jac->use_true_matrix = PETSC_FALSE;
275: jac->its = 0;
277: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCKSPSetUseTrue_C","PCKSPSetUseTrue_KSP",
278: PCKSPSetUseTrue_KSP);
279: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCKSPGetKSP_C","PCKSPGetKSP_KSP",
280: PCKSPGetKSP_KSP);
282: return(0);
283: }