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;
15: static PetscErrorCode PCKSPCreateKSP_KSP(PC pc)
16: {
18: const char *prefix;
19: PC_KSP *jac = (PC_KSP *)pc->data;
22: KSPCreate(((PetscObject)pc)->comm,&jac->ksp);
23: PetscObjectIncrementTabLevel((PetscObject)jac->ksp,(PetscObject)pc,1);
24: PCGetOptionsPrefix(pc,&prefix);
25: KSPSetOptionsPrefix(jac->ksp,prefix);
26: KSPAppendOptionsPrefix(jac->ksp,"ksp_");
27: return(0);
28: }
32: static PetscErrorCode PCApply_KSP(PC pc,Vec x,Vec y)
33: {
35: PetscInt its;
36: PC_KSP *jac = (PC_KSP*)pc->data;
39: KSPSetInitialGuessNonzero(jac->ksp,pc->nonzero_guess);
40: KSPSolve(jac->ksp,x,y);
41: KSPGetIterationNumber(jac->ksp,&its);
42: jac->its += its;
43: return(0);
44: }
48: static PetscErrorCode PCApplyTranspose_KSP(PC pc,Vec x,Vec y)
49: {
51: PetscInt its;
52: PC_KSP *jac = (PC_KSP*)pc->data;
55: KSPSolveTranspose(jac->ksp,x,y);
56: KSPGetIterationNumber(jac->ksp,&its);
57: jac->its += its;
58: return(0);
59: }
63: static PetscErrorCode PCSetUp_KSP(PC pc)
64: {
66: PC_KSP *jac = (PC_KSP*)pc->data;
67: Mat mat;
68: PetscTruth A;
71: if (!jac->ksp) {PCKSPCreateKSP_KSP(pc);}
72: KSPSetFromOptions(jac->ksp);
73: if (jac->use_true_matrix) mat = pc->mat;
74: else mat = pc->pmat;
76: KSPGetOperatorsSet(jac->ksp,&A,PETSC_NULL);
77: if (!A) {
78: KSPSetOperators(jac->ksp,mat,pc->pmat,pc->flag);
79: }
80: KSPSetUp(jac->ksp);
81: return(0);
82: }
84: /* Default destroy, if it has never been setup */
87: static PetscErrorCode PCDestroy_KSP(PC pc)
88: {
89: PC_KSP *jac = (PC_KSP*)pc->data;
93: if (jac->ksp) {KSPDestroy(jac->ksp);}
94: PetscFree(jac);
95: return(0);
96: }
100: static PetscErrorCode PCView_KSP(PC pc,PetscViewer viewer)
101: {
102: PC_KSP *jac = (PC_KSP*)pc->data;
104: PetscTruth iascii;
107: if (!jac->ksp) {PCKSPCreateKSP_KSP(pc);}
108: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
109: if (iascii) {
110: if (jac->use_true_matrix) {
111: PetscViewerASCIIPrintf(viewer,"Using true matrix (not preconditioner matrix) on inner solve\n");
112: }
113: PetscViewerASCIIPrintf(viewer,"KSP and PC on KSP preconditioner follow\n");
114: PetscViewerASCIIPrintf(viewer,"---------------------------------\n");
115: } else {
116: SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for this object",((PetscObject)viewer)->type_name);
117: }
118: PetscViewerASCIIPushTab(viewer);
119: KSPView(jac->ksp,viewer);
120: PetscViewerASCIIPopTab(viewer);
121: if (iascii) {
122: PetscViewerASCIIPrintf(viewer,"---------------------------------\n");
123: }
124: return(0);
125: }
129: static PetscErrorCode PCSetFromOptions_KSP(PC pc)
130: {
132: PetscTruth flg;
135: PetscOptionsHead("KSP preconditioner options");
136: PetscOptionsName("-pc_ksp_true","Use true matrix to define inner linear system, not preconditioner matrix","PCKSPSetUseTrue",&flg);
137: if (flg) {
138: PCKSPSetUseTrue(pc);
139: }
140: PetscOptionsTail();
141: return(0);
142: }
144: /* ----------------------------------------------------------------------------------*/
149: PetscErrorCode PCKSPSetUseTrue_KSP(PC pc)
150: {
151: PC_KSP *jac;
154: jac = (PC_KSP*)pc->data;
155: jac->use_true_matrix = PETSC_TRUE;
156: return(0);
157: }
163: PetscErrorCode PCKSPGetKSP_KSP(PC pc,KSP *ksp)
164: {
165: PC_KSP *jac = (PC_KSP*)pc->data;
169: if (!jac->ksp) {PCKSPCreateKSP_KSP(pc);}
170: *ksp = jac->ksp;
171: return(0);
172: }
177: /*@
178: PCKSPSetUseTrue - Sets a flag to indicate that the true matrix (rather than
179: the matrix used to define the preconditioner) is used to compute the
180: residual inside the inner solve.
182: Collective on PC
184: Input Parameters:
185: . pc - the preconditioner context
187: Options Database Key:
188: . -pc_ksp_true - Activates PCKSPSetUseTrue()
190: Note:
191: For the common case in which the preconditioning and linear
192: system matrices are identical, this routine is unnecessary.
194: Level: advanced
196: .keywords: PC, KSP, set, true, local, flag
198: .seealso: PCSetOperators(), PCBJacobiSetUseTrueLocal()
199: @*/
200: PetscErrorCode PCKSPSetUseTrue(PC pc)
201: {
202: PetscErrorCode ierr,(*f)(PC);
206: PetscObjectQueryFunction((PetscObject)pc,"PCKSPSetUseTrue_C",(void (**)(void))&f);
207: if (f) {
208: (*f)(pc);
209: }
210: return(0);
211: }
215: /*@
216: PCKSPGetKSP - Gets the KSP context for a KSP PC.
218: Not Collective but KSP returned is parallel if PC was parallel
220: Input Parameter:
221: . pc - the preconditioner context
223: Output Parameters:
224: . ksp - the PC solver
226: Notes:
227: You must call KSPSetUp() before calling PCKSPGetKSP().
229: Level: advanced
231: .keywords: PC, KSP, get, context
232: @*/
233: PetscErrorCode PCKSPGetKSP(PC pc,KSP *ksp)
234: {
235: PetscErrorCode ierr,(*f)(PC,KSP*);
240: PetscObjectQueryFunction((PetscObject)pc,"PCKSPGetKSP_C",(void (**)(void))&f);
241: if (f) {
242: (*f)(pc,ksp);
243: }
244: return(0);
245: }
247: /* ----------------------------------------------------------------------------------*/
249: /*MC
250: PCKSP - Defines a preconditioner that can consist of any KSP solver.
251: This allows, for example, embedding a Krylov method inside a preconditioner.
253: Options Database Key:
254: . -pc_ksp_true - use the matrix that defines the linear system as the matrix for the
255: inner solver, otherwise by default it uses the matrix used to construct
256: the preconditioner (see PCSetOperators())
258: Level: intermediate
260: Concepts: inner iteration
262: Notes: Using a Krylov method inside another Krylov method can be dangerous (you get divergence or
263: the incorrect answer) unless you use KSPFGMRES as the other Krylov method
266: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
267: PCSHELL, PCCOMPOSITE, PCKSPUseTrue(), PCKSPGetKSP()
269: M*/
274: PetscErrorCode PCCreate_KSP(PC pc)
275: {
277: PC_KSP *jac;
280: PetscNewLog(pc,PC_KSP,&jac);
281: pc->ops->apply = PCApply_KSP;
282: pc->ops->applytranspose = PCApplyTranspose_KSP;
283: pc->ops->setup = PCSetUp_KSP;
284: pc->ops->destroy = PCDestroy_KSP;
285: pc->ops->setfromoptions = PCSetFromOptions_KSP;
286: pc->ops->view = PCView_KSP;
287: pc->ops->applyrichardson = 0;
289: pc->data = (void*)jac;
290:
292: jac->use_true_matrix = PETSC_FALSE;
293: jac->its = 0;
295: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCKSPSetUseTrue_C","PCKSPSetUseTrue_KSP",
296: PCKSPSetUseTrue_KSP);
297: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCKSPGetKSP_C","PCKSPGetKSP_KSP",
298: PCKSPGetKSP_KSP);
300: return(0);
301: }