Actual source code: pcksp.c
1: /*$Id: pcsles.c,v 1.39 2001/04/10 19:36:17 bsmith Exp $*/
4: #include src/ksp/pc/pcimpl.h
5: #include petscksp.h
7: typedef struct {
8: PetscTruth use_true_matrix; /* use mat rather than pmat in inner linear solve */
9: KSP ksp;
10: int its; /* total number of iterations KSP uses */
11: } PC_KSP;
15: static int PCApply_KSP(PC pc,Vec x,Vec y)
16: {
17: int ierr,its;
18: PC_KSP *jac = (PC_KSP*)pc->data;
21: KSPSetRhs(jac->ksp,x);
22: KSPSetSolution(jac->ksp,y);
23: KSPSolve(jac->ksp);
24: KSPGetIterationNumber(jac->ksp,&its);
25: jac->its += its;
26: return(0);
27: }
31: static int PCApplyTranspose_KSP(PC pc,Vec x,Vec y)
32: {
33: int its,ierr;
34: PC_KSP *jac = (PC_KSP*)pc->data;
37: KSPSetRhs(jac->ksp,x);
38: KSPSetSolution(jac->ksp,y);
39: KSPSolveTranspose(jac->ksp);
40: KSPGetIterationNumber(jac->ksp,&its);
41: jac->its += its;
42: return(0);
43: }
47: static int PCSetUp_KSP(PC pc)
48: {
49: int ierr;
50: PC_KSP *jac = (PC_KSP*)pc->data;
51: Mat mat;
54: KSPSetFromOptions(jac->ksp);
55: if (jac->use_true_matrix) mat = pc->mat;
56: else mat = pc->pmat;
58: KSPSetOperators(jac->ksp,mat,pc->pmat,pc->flag);
59: KSPSetRhs(jac->ksp,pc->vec);
60: KSPSetSolution(jac->ksp,pc->vec);
61: KSPSetUp(jac->ksp);
62: return(0);
63: }
65: /* Default destroy, if it has never been setup */
68: static int PCDestroy_KSP(PC pc)
69: {
70: PC_KSP *jac = (PC_KSP*)pc->data;
71: int ierr;
74: KSPDestroy(jac->ksp);
75: PetscFree(jac);
76: return(0);
77: }
81: static int PCView_KSP(PC pc,PetscViewer viewer)
82: {
83: PC_KSP *jac = (PC_KSP*)pc->data;
84: int ierr;
85: PetscTruth isascii;
88: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
89: if (isascii) {
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(1,"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 (isascii) {
102: PetscViewerASCIIPrintf(viewer,"---------------------------------\n");
103: }
104: return(0);
105: }
109: static int PCSetFromOptions_KSP(PC pc){
110: int ierr;
111: PetscTruth flg;
114: PetscOptionsHead("KSP preconditioner options");
115: PetscOptionsName("-pc_ksp_true","Use true matrix to define inner linear system, not preconditioner matrix","PCKSPSetUseTrue",&flg);
116: if (flg) {
117: PCKSPSetUseTrue(pc);
118: }
119: PetscOptionsTail();
120: return(0);
121: }
123: /* ----------------------------------------------------------------------------------*/
125: EXTERN_C_BEGIN
128: int PCKSPSetUseTrue_KSP(PC pc)
129: {
130: PC_KSP *jac;
133: jac = (PC_KSP*)pc->data;
134: jac->use_true_matrix = PETSC_TRUE;
135: return(0);
136: }
137: EXTERN_C_END
139: EXTERN_C_BEGIN
142: int PCKSPGetKSP_KSP(PC pc,KSP *ksp)
143: {
144: PC_KSP *jac;
147: jac = (PC_KSP*)pc->data;
148: *ksp = jac->ksp;
149: return(0);
150: }
151: EXTERN_C_END
155: /*@
156: PCKSPSetUseTrue - Sets a flag to indicate that the true matrix (rather than
157: the matrix used to define the preconditioner) is used to compute the
158: residual inside the inner solve.
160: Collective on PC
162: Input Parameters:
163: . pc - the preconditioner context
165: Options Database Key:
166: . -pc_ksp_true - Activates PCKSPSetUseTrue()
168: Note:
169: For the common case in which the preconditioning and linear
170: system matrices are identical, this routine is unnecessary.
172: Level: advanced
174: .keywords: PC, KSP, set, true, local, flag
176: .seealso: PCSetOperators(), PCBJacobiSetUseTrueLocal()
177: @*/
178: int PCKSPSetUseTrue(PC pc)
179: {
180: int ierr,(*f)(PC);
184: PetscObjectQueryFunction((PetscObject)pc,"PCKSPSetUseTrue_C",(void (**)(void))&f);
185: if (f) {
186: (*f)(pc);
187: }
188: return(0);
189: }
193: /*@C
194: PCKSPGetKSP - Gets the KSP context for a KSP PC.
196: Not Collective but KSP returned is parallel if PC was parallel
198: Input Parameter:
199: . pc - the preconditioner context
201: Output Parameters:
202: . ksp - the PC solver
204: Notes:
205: You must call KSPSetUp() before calling PCKSPGetKSP().
207: Level: advanced
209: .keywords: PC, KSP, get, context
210: @*/
211: int PCKSPGetKSP(PC pc,KSP *ksp)
212: {
213: int ierr,(*f)(PC,KSP*);
218: if (!pc->setupcalled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call KSPSetUp first");
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*/
250: EXTERN_C_BEGIN
253: int PCCreate_KSP(PC pc)
254: {
255: int ierr;
256: char *prefix;
257: PC_KSP *jac;
260: PetscNew(PC_KSP,&jac);
261: PetscLogObjectMemory(pc,sizeof(PC_KSP));
262: pc->ops->apply = PCApply_KSP;
263: pc->ops->applytranspose = PCApplyTranspose_KSP;
264: pc->ops->setup = PCSetUp_KSP;
265: pc->ops->destroy = PCDestroy_KSP;
266: pc->ops->setfromoptions = PCSetFromOptions_KSP;
267: pc->ops->view = PCView_KSP;
268: pc->ops->applyrichardson = 0;
270: pc->data = (void*)jac;
271: KSPCreate(pc->comm,&jac->ksp);
273: PCGetOptionsPrefix(pc,&prefix);
274: KSPSetOptionsPrefix(jac->ksp,prefix);
275: KSPAppendOptionsPrefix(jac->ksp,"ksp_");
276: jac->use_true_matrix = PETSC_FALSE;
277: jac->its = 0;
279: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCKSPSetUseTrue_C","PCKSPSetUseTrue_KSP",
280: PCKSPSetUseTrue_KSP);
281: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCKSPGetKSP_C","PCKSPGetKSP_KSP",
282: PCKSPGetKSP_KSP);
284: return(0);
285: }
286: EXTERN_C_END