xref: /petsc/src/ksp/ksp/impls/fetidp/fetidp.c (revision 43b14c2c43eccf797c1bca6e0ea443d70841cdc8)
1 #include <petsc/private/kspimpl.h> /*I <petscksp.h> I*/
2 #include <../src/ksp/pc/impls/bddc/bddc.h>
3 
4 /*
5     This file implements the FETI-DP method in PETSc as part of KSP.
6 */
7 typedef struct {
8   KSP parentksp;
9 } KSP_FETIDPMon;
10 
11 typedef struct {
12   KSP           innerksp;        /* the KSP for the Lagrange multipliers */
13   PC            innerbddc;       /* the inner BDDC object */
14   PetscBool     fully_redundant; /* true for using a fully redundant set of multipliers */
15   PetscBool     userbddc;        /* true if the user provided the PCBDDC object */
16   KSP_FETIDPMon *monctx;         /* monitor context, used to pass user defined monitors for the physical */
17 } KSP_FETIDP;
18 
19 #undef __FUNCT__
20 #define __FUNCT__ "KSPFETIDPGetInnerKSP_FETIDP"
21 static PetscErrorCode KSPFETIDPGetInnerKSP_FETIDP(KSP ksp,KSP* innerksp)
22 {
23   KSP_FETIDP     *fetidp = (KSP_FETIDP*)ksp->data;
24 
25   PetscFunctionBegin;
26   *innerksp = fetidp->innerksp;
27   PetscFunctionReturn(0);
28 }
29 
30 #undef __FUNCT__
31 #define __FUNCT__ "KSPFETIDPGetInnerKSP"
32 /*@
33  KSPFETIDPGetInnerKSP - Get the KSP object for the Lagrange multipliers
34 
35    Input Parameters:
36 +  ksp - the FETI-DP KSP
37 -  innerksp - the KSP for the multipliers
38 
39    Level: advanced
40 
41    Notes:
42 
43 .seealso: MATIS, PCBDDC, KSPFETIDPSetInnerBDDC, KSPFETIDPGetInnerBDDC
44 @*/
45 PetscErrorCode KSPFETIDPGetInnerKSP(KSP ksp,KSP* innerksp)
46 {
47   PetscErrorCode ierr;
48 
49   PetscFunctionBegin;
50   PetscValidHeaderSpecific(ksp,KSP_CLASSID,1);
51   PetscValidPointer(innerksp,2);
52   ierr = PetscUseMethod(ksp,"KSPFETIDPGetInnerKSP_C",(KSP,KSP*),(ksp,innerksp));CHKERRQ(ierr);
53   PetscFunctionReturn(0);
54 }
55 
56 #undef __FUNCT__
57 #define __FUNCT__ "KSPFETIDPGetInnerBDDC_FETIDP"
58 static PetscErrorCode KSPFETIDPGetInnerBDDC_FETIDP(KSP ksp,PC* pc)
59 {
60   KSP_FETIDP     *fetidp = (KSP_FETIDP*)ksp->data;
61 
62   PetscFunctionBegin;
63   *pc = fetidp->innerbddc;
64   PetscFunctionReturn(0);
65 }
66 
67 #undef __FUNCT__
68 #define __FUNCT__ "KSPFETIDPGetInnerBDDC"
69 /*@
70  KSPFETIDPGetInnerBDDC - Get the PCBDDC object used to setup the FETI-DP matrix for the Lagrange multipliers
71 
72    Input Parameters:
73 +  ksp - the FETI-DP KSP
74 -  pc - the PCBDDC object
75 
76    Level: advanced
77 
78    Notes:
79 
80 .seealso: MATIS, PCBDDC, KSPFETIDPSetInnerBDDC, KSPFETIDPGetInnerKSP
81 @*/
82 PetscErrorCode KSPFETIDPGetInnerBDDC(KSP ksp,PC* pc)
83 {
84   PetscErrorCode ierr;
85 
86   PetscFunctionBegin;
87   PetscValidHeaderSpecific(ksp,KSP_CLASSID,1);
88   PetscValidPointer(pc,2);
89   ierr = PetscUseMethod(ksp,"KSPFETIDPGetInnerBDDC_C",(KSP,PC*),(ksp,pc));CHKERRQ(ierr);
90   PetscFunctionReturn(0);
91 }
92 
93 #undef __FUNCT__
94 #define __FUNCT__ "KSPFETIDPSetInnerBDDC_FETIDP"
95 static PetscErrorCode KSPFETIDPSetInnerBDDC_FETIDP(KSP ksp,PC pc)
96 {
97   KSP_FETIDP     *fetidp = (KSP_FETIDP*)ksp->data;
98   PetscErrorCode ierr;
99 
100   PetscFunctionBegin;
101   ierr = PetscObjectReference((PetscObject)pc);CHKERRQ(ierr);
102   ierr = PCDestroy(&fetidp->innerbddc);CHKERRQ(ierr);
103   fetidp->innerbddc = pc;
104   fetidp->userbddc = PETSC_TRUE;
105   PetscFunctionReturn(0);
106 }
107 
108 #undef __FUNCT__
109 #define __FUNCT__ "KSPFETIDPSetInnerBDDC"
110 /*@
111  KSPFETIDPSetInnerBDDC - Set the PCBDDC object used to setup the FETI-DP matrix for the Lagrange multipliers
112 
113    Collective on KSP
114 
115    Input Parameters:
116 +  ksp - the FETI-DP KSP
117 -  pc - the PCBDDC object
118 
119    Level: advanced
120 
121    Notes:
122 
123 .seealso: MATIS, PCBDDC, KSPFETIDPGetInnerBDDC, KSPFETIDPGetInnerKSP
124 @*/
125 PetscErrorCode KSPFETIDPSetInnerBDDC(KSP ksp,PC pc)
126 {
127   PetscBool      isbddc;
128   PetscErrorCode ierr;
129 
130   PetscFunctionBegin;
131   PetscValidHeaderSpecific(ksp,KSP_CLASSID,1);
132   PetscValidHeaderSpecific(pc,PC_CLASSID,2);
133   ierr = PetscObjectTypeCompare((PetscObject)pc,PCBDDC,&isbddc);CHKERRQ(ierr);
134   if (!isbddc) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONG,"KSPFETIDPSetInnerBDDC need a PCBDDC preconditioner");
135   ierr = PetscTryMethod(ksp,"KSPFETIDPSetInnerBDDC_C",(KSP,PC),(ksp,pc));CHKERRQ(ierr);
136   PetscFunctionReturn(0);
137 }
138 
139 #undef __FUNCT__
140 #define __FUNCT__ "KSPBuildSolution_FETIDP"
141 static PetscErrorCode KSPBuildSolution_FETIDP(KSP ksp,Vec v,Vec *V)
142 {
143   KSP_FETIDP     *fetidp = (KSP_FETIDP*)ksp->data;
144   Mat            F;
145   Vec            Xl;
146   PetscErrorCode ierr;
147 
148   PetscFunctionBegin;
149   ierr = KSPGetOperators(fetidp->innerksp,&F,NULL);CHKERRQ(ierr);
150   ierr = KSPBuildSolution(fetidp->innerksp,NULL,&Xl);CHKERRQ(ierr);
151   if (v) {
152     ierr = PCBDDCMatFETIDPGetSolution(F,Xl,v);CHKERRQ(ierr);
153     *V = v;
154   } else {
155     ierr = PCBDDCMatFETIDPGetSolution(F,Xl,*V);CHKERRQ(ierr);
156   }
157   PetscFunctionReturn(0);
158 }
159 
160 #undef __FUNCT__
161 #define __FUNCT__ "KSPMonitor_FETIDP"
162 static PetscErrorCode KSPMonitor_FETIDP(KSP ksp,PetscInt it,PetscReal rnorm,void* ctx)
163 {
164   KSP_FETIDPMon  *monctx = (KSP_FETIDPMon*)ctx;
165   PetscErrorCode ierr;
166 
167   PetscFunctionBegin;
168   ierr = KSPMonitor(monctx->parentksp,it,rnorm);CHKERRQ(ierr);
169   PetscFunctionReturn(0);
170 }
171 
172 #undef __FUNCT__
173 #define __FUNCT__ "KSPComputeEigenvalues_FETIDP"
174 static PetscErrorCode KSPComputeEigenvalues_FETIDP(KSP ksp,PetscInt nmax,PetscReal *r,PetscReal *c,PetscInt *neig)
175 {
176   KSP_FETIDP     *fetidp = (KSP_FETIDP*)ksp->data;
177   PetscErrorCode ierr;
178 
179   PetscFunctionBegin;
180   ierr = KSPComputeEigenvalues(fetidp->innerksp,nmax,r,c,neig);CHKERRQ(ierr);
181   PetscFunctionReturn(0);
182 }
183 
184 #undef __FUNCT__
185 #define __FUNCT__ "KSPComputeExtremeSingularValues_FETIDP"
186 static PetscErrorCode KSPComputeExtremeSingularValues_FETIDP(KSP ksp,PetscReal *emax,PetscReal *emin)
187 {
188   KSP_FETIDP     *fetidp = (KSP_FETIDP*)ksp->data;
189   PetscErrorCode ierr;
190 
191   PetscFunctionBegin;
192   ierr = KSPComputeExtremeSingularValues(fetidp->innerksp,emax,emin);CHKERRQ(ierr);
193   PetscFunctionReturn(0);
194 }
195 
196 #undef __FUNCT__
197 #define __FUNCT__ "KSPSetUp_FETIDP"
198 static PetscErrorCode KSPSetUp_FETIDP(KSP ksp)
199 {
200   KSP_FETIDP     *fetidp = (KSP_FETIDP*)ksp->data;
201   PC_BDDC        *pcbddc = (PC_BDDC*)fetidp->innerbddc->data;
202   Mat            A,Ap;
203   PetscBool      flg;
204   PetscErrorCode ierr;
205 
206   PetscFunctionBegin;
207   /* set up BDDC */
208   ierr = KSPGetOperators(ksp,&A,&Ap);CHKERRQ(ierr);
209   ierr = PCSetOperators(fetidp->innerbddc,A,Ap);CHKERRQ(ierr);
210   ierr = PCSetUp(fetidp->innerbddc);CHKERRQ(ierr);
211 
212   /* FETI-DP as it is implemented needs an exact coarse solver */
213   if (pcbddc->coarse_ksp) {
214     ierr = KSPSetTolerances(pcbddc->coarse_ksp,PETSC_SMALL,PETSC_SMALL,PETSC_DEFAULT,1000);CHKERRQ(ierr);
215     ierr = KSPSetNormType(pcbddc->coarse_ksp,KSP_NORM_DEFAULT);CHKERRQ(ierr);
216   }
217   /* FETI-DP as it is implemented needs exact local solvers */
218   ierr = KSPSetTolerances(pcbddc->ksp_R,PETSC_SMALL,PETSC_SMALL,PETSC_DEFAULT,1000);CHKERRQ(ierr);
219   ierr = KSPSetNormType(pcbddc->ksp_R,KSP_NORM_DEFAULT);CHKERRQ(ierr);
220   /* if the primal space is changed, setup F */
221   if (pcbddc->new_primal_space) {
222     Mat F; /* the FETI-DP matrix */
223     PC  D; /* the FETI-DP preconditioner */
224     ierr = KSPReset(fetidp->innerksp);CHKERRQ(ierr);
225     ierr = PCBDDCCreateFETIDPOperators(fetidp->innerbddc,fetidp->fully_redundant,&F,&D);CHKERRQ(ierr);
226     ierr = KSPSetOperators(fetidp->innerksp,F,F);CHKERRQ(ierr);
227     ierr = KSPSetTolerances(fetidp->innerksp,ksp->rtol,ksp->abstol,ksp->divtol,ksp->max_it);CHKERRQ(ierr);
228     ierr = KSPSetPC(fetidp->innerksp,D);CHKERRQ(ierr);
229     ierr = MatCreateVecs(F,&(fetidp->innerksp)->vec_rhs,&(fetidp->innerksp)->vec_sol);CHKERRQ(ierr);
230     ierr = MatDestroy(&F);CHKERRQ(ierr);
231     ierr = PCDestroy(&D);CHKERRQ(ierr);
232   }
233   /* propagate settings to inner solve */
234   ierr = KSPGetComputeSingularValues(ksp,&flg);CHKERRQ(ierr);
235   ierr = KSPSetComputeSingularValues(fetidp->innerksp,flg);CHKERRQ(ierr);
236   if (ksp->res_hist) {
237     ierr = KSPSetResidualHistory(fetidp->innerksp,ksp->res_hist,ksp->res_hist_max,ksp->res_hist_reset);CHKERRQ(ierr);
238   }
239   ierr = KSPSetUp(fetidp->innerksp);CHKERRQ(ierr);
240   PetscFunctionReturn(0);
241 }
242 
243 #undef __FUNCT__
244 #define __FUNCT__ "KSPSolve_FETIDP"
245 static PetscErrorCode KSPSolve_FETIDP(KSP ksp)
246 {
247   PetscErrorCode ierr;
248   Mat            F;
249   Vec            X,B,Xl,Bl;
250   KSP_FETIDP     *fetidp = (KSP_FETIDP*)ksp->data;
251 
252   PetscFunctionBegin;
253   ierr = KSPGetRhs(ksp,&B);CHKERRQ(ierr);
254   ierr = KSPGetSolution(ksp,&X);CHKERRQ(ierr);
255   ierr = KSPGetOperators(fetidp->innerksp,&F,NULL);CHKERRQ(ierr);
256   ierr = KSPGetRhs(fetidp->innerksp,&Bl);CHKERRQ(ierr);
257   ierr = KSPGetSolution(fetidp->innerksp,&Xl);CHKERRQ(ierr);
258   ierr = PCBDDCMatFETIDPGetRHS(F,B,Bl);CHKERRQ(ierr);
259   if (ksp->transpose_solve) {
260     ierr = KSPSolveTranspose(fetidp->innerksp,Bl,Xl);CHKERRQ(ierr);
261   } else {
262     ierr = KSPSolve(fetidp->innerksp,Bl,Xl);CHKERRQ(ierr);
263   }
264   ierr = PCBDDCMatFETIDPGetSolution(F,Xl,X);CHKERRQ(ierr);
265   /* update ksp with stats from inner ksp */
266   ierr = KSPGetConvergedReason(fetidp->innerksp,&ksp->reason);CHKERRQ(ierr);
267   ierr = KSPGetIterationNumber(fetidp->innerksp,&ksp->its);CHKERRQ(ierr);
268   ksp->totalits += ksp->its;
269   ierr = KSPGetResidualHistory(fetidp->innerksp,NULL,&ksp->res_hist_len);CHKERRQ(ierr);
270   PetscFunctionReturn(0);
271 }
272 
273 #undef __FUNCT__
274 #define __FUNCT__ "KSPDestroy_FETIDP"
275 static PetscErrorCode KSPDestroy_FETIDP(KSP ksp)
276 {
277   KSP_FETIDP     *fetidp = (KSP_FETIDP*)ksp->data;
278   PetscErrorCode ierr;
279 
280   PetscFunctionBegin;
281   ierr = PCDestroy(&fetidp->innerbddc);CHKERRQ(ierr);
282   ierr = KSPDestroy(&fetidp->innerksp);CHKERRQ(ierr);
283   ierr = PetscFree(fetidp->monctx);CHKERRQ(ierr);
284   ierr = PetscObjectComposeFunction((PetscObject)ksp,"KSPFETIDPSetInnerBDDC_C",NULL);CHKERRQ(ierr);
285   ierr = PetscObjectComposeFunction((PetscObject)ksp,"KSPFETIDPGetInnerBDDC_C",NULL);CHKERRQ(ierr);
286   ierr = PetscObjectComposeFunction((PetscObject)ksp,"KSPFETIDPGetInnerKSP_C",NULL);CHKERRQ(ierr);
287   ierr = PetscFree(ksp->data);CHKERRQ(ierr);
288   PetscFunctionReturn(0);
289 }
290 
291 #undef __FUNCT__
292 #define __FUNCT__ "KSPView_FETIDP"
293 static PetscErrorCode KSPView_FETIDP(KSP ksp,PetscViewer viewer)
294 {
295   KSP_FETIDP     *fetidp = (KSP_FETIDP*)ksp->data;
296   PetscErrorCode ierr;
297   PetscBool      iascii;
298 
299   PetscFunctionBegin;
300   ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr);
301   if (iascii) {
302     ierr = PetscViewerASCIIPrintf(viewer,"  FETI_DP: fully redundant: %d\n",fetidp->fully_redundant);CHKERRQ(ierr);
303     ierr = PetscViewerASCIIPrintf(viewer,"  FETI_DP: inner solver details\n");CHKERRQ(ierr);
304     ierr = PetscViewerASCIIAddTab(viewer,2);CHKERRQ(ierr);
305   }
306   ierr = KSPView(fetidp->innerksp,viewer);CHKERRQ(ierr);
307   if (iascii) {
308     ierr = PetscViewerASCIISubtractTab(viewer,2);CHKERRQ(ierr);
309     ierr = PetscViewerASCIIPrintf(viewer,"  FETI_DP: BDDC solver details\n");CHKERRQ(ierr);
310     ierr = PetscViewerASCIIAddTab(viewer,2);CHKERRQ(ierr);
311   }
312   ierr = PCView(fetidp->innerbddc,viewer);CHKERRQ(ierr);
313   if (iascii) {
314     ierr = PetscViewerASCIISubtractTab(viewer,2);CHKERRQ(ierr);
315   }
316   PetscFunctionReturn(0);
317 }
318 
319 #undef __FUNCT__
320 #define __FUNCT__ "KSPSetFromOptions_FETIDP"
321 static PetscErrorCode KSPSetFromOptions_FETIDP(PetscOptionItems *PetscOptionsObject,KSP ksp)
322 {
323   KSP_FETIDP     *fetidp = (KSP_FETIDP*)ksp->data;
324   PetscErrorCode ierr;
325 
326   PetscFunctionBegin;
327   /* set options prefixes for the inner objects, since the parent prefix will be valid at this point */
328   ierr = PetscObjectSetOptionsPrefix((PetscObject)fetidp->innerksp,((PetscObject)ksp)->prefix);CHKERRQ(ierr);
329   ierr = PetscObjectAppendOptionsPrefix((PetscObject)fetidp->innerksp,"fetidp_");CHKERRQ(ierr);
330   if (!fetidp->userbddc) {
331     ierr = PetscObjectSetOptionsPrefix((PetscObject)fetidp->innerbddc,((PetscObject)ksp)->prefix);CHKERRQ(ierr);
332     ierr = PetscObjectAppendOptionsPrefix((PetscObject)fetidp->innerbddc,"fetidp_bddc_");CHKERRQ(ierr);
333   }
334   ierr = PetscOptionsHead(PetscOptionsObject,"KSP FETIDP options");CHKERRQ(ierr);
335   ierr = PetscOptionsBool("-ksp_fetidp_fullyredundant","Use fully redundant multipliers","none",fetidp->fully_redundant,&fetidp->fully_redundant,NULL);CHKERRQ(ierr);
336   ierr = PetscOptionsTail();CHKERRQ(ierr);
337   ierr = PCSetFromOptions(fetidp->innerbddc);CHKERRQ(ierr);
338   ierr = KSPSetFromOptions(fetidp->innerksp);CHKERRQ(ierr);
339   PetscFunctionReturn(0);
340 }
341 
342 /*MC
343      KSPFETIDP - The FETI-DP method
344 
345    This class implements the FETI-DP method [1].
346    The preconditioning matrix for the KSP must be of type MATIS.
347    The FETI-DP linear system (automatically generated constructing an internal PCBDDC object) is solved using an inner KSP object.
348 
349    Options Database Keys:
350 .   -ksp_fetidp_fullyredundant <false> : use a fully redundant set of Lagrange multipliers
351 
352    Level: Advanced
353 
354    Notes: Options for the inner KSP and for the customization of the PCBDDC object can be specified at command line by using the prefixes -fetidp_ and -fetidp_bddc_. E.g.,
355 .vb
356       -fetidp_ksp_type gmres -fetidp_bddc_pc_bddc_symmetric false
357 .ve
358    will use GMRES for the solution of the linear system on the Lagrange multipliers, generated using a non-symmetric PCBDDC.
359 
360    References:
361 .  [1] - C. Farhat, M. Lesoinne, P. LeTallec, K. Pierson, and D. Rixen, FETI-DP: a dual-primal unified FETI method. I. A faster alternative to the two-level FETI method, Internat. J. Numer. Methods Engrg., 50 (2001), pp. 1523--1544
362 
363 .seealso: MATIS, PCBDDC, KSPFETIDPSetInnerBDDC, KSPFETIDPGetInnerBDDC, KSPFETIDPGetInnerKSP
364 M*/
365 #undef __FUNCT__
366 #define __FUNCT__ "KSPCreate_FETIDP"
367 PETSC_EXTERN PetscErrorCode KSPCreate_FETIDP(KSP ksp)
368 {
369   PetscErrorCode ierr;
370   KSP_FETIDP     *fetidp;
371   KSP_FETIDPMon  *monctx;
372   PC_BDDC        *pcbddc;
373   PC             pc;
374 
375   PetscFunctionBegin;
376   ierr = PetscNewLog(ksp,&fetidp);CHKERRQ(ierr);
377   ksp->data = (void*)fetidp;
378   ksp->ops->setup                        = KSPSetUp_FETIDP;
379   ksp->ops->solve                        = KSPSolve_FETIDP;
380   ksp->ops->destroy                      = KSPDestroy_FETIDP;
381   ksp->ops->computeeigenvalues           = KSPComputeEigenvalues_FETIDP;
382   ksp->ops->computeextremesingularvalues = KSPComputeExtremeSingularValues_FETIDP;
383   ksp->ops->view                         = KSPView_FETIDP;
384   ksp->ops->setfromoptions               = KSPSetFromOptions_FETIDP;
385   ksp->ops->buildsolution                = KSPBuildSolution_FETIDP;
386   ksp->ops->buildresidual                = KSPBuildResidualDefault;
387   /* create the inner KSP for the Lagrange multipliers */
388   ierr = KSPCreate(PetscObjectComm((PetscObject)ksp),&fetidp->innerksp);CHKERRQ(ierr);
389   ierr = KSPGetPC(fetidp->innerksp,&pc);CHKERRQ(ierr);
390   ierr = PCSetType(pc,PCNONE);CHKERRQ(ierr);
391   ierr = PetscLogObjectParent((PetscObject)ksp,(PetscObject)fetidp->innerksp);CHKERRQ(ierr);
392   /* monitor */
393   ierr = PetscNew(&monctx);CHKERRQ(ierr);
394   monctx->parentksp = ksp;
395   fetidp->monctx = monctx;
396   ierr = KSPMonitorSet(fetidp->innerksp,KSPMonitor_FETIDP,fetidp->monctx,NULL);CHKERRQ(ierr);
397   /* create the inner BDDC */
398   ierr = PCCreate(PetscObjectComm((PetscObject)ksp),&fetidp->innerbddc);CHKERRQ(ierr);
399   ierr = PCSetType(fetidp->innerbddc,PCBDDC);CHKERRQ(ierr);
400   /* make sure we always obtain a consistent FETI-DP matrix
401      for symmetric problems, the user can always customize it through the command line */
402   pcbddc = (PC_BDDC*)fetidp->innerbddc->data;
403   pcbddc->symmetric_primal = PETSC_FALSE;
404   ierr = PetscLogObjectParent((PetscObject)ksp,(PetscObject)fetidp->innerbddc);CHKERRQ(ierr);
405   /* composed functions */
406   ierr = PetscObjectComposeFunction((PetscObject)ksp,"KSPFETIDPSetInnerBDDC_C",KSPFETIDPSetInnerBDDC_FETIDP);CHKERRQ(ierr);
407   ierr = PetscObjectComposeFunction((PetscObject)ksp,"KSPFETIDPGetInnerBDDC_C",KSPFETIDPGetInnerBDDC_FETIDP);CHKERRQ(ierr);
408   ierr = PetscObjectComposeFunction((PetscObject)ksp,"KSPFETIDPGetInnerKSP_C",KSPFETIDPGetInnerKSP_FETIDP);CHKERRQ(ierr);
409   /* need to call KSPSetUp_FETIDP even with KSP_SETUP_NEWMATRIX */
410   ksp->setupnewmatrix = PETSC_TRUE;
411   PetscFunctionReturn(0);
412 }
413