xref: /petsc/src/ksp/pc/impls/factor/cholesky/cholesky.c (revision 3fa76a5b188a02dd677974deca162194da798f66)
1 #define PETSCKSP_DLL
2 
3 /*
4    Defines a direct factorization preconditioner for any Mat implementation
5    Note: this need not be consided a preconditioner since it supplies
6          a direct solver.
7 */
8 #include "private/pcimpl.h"                /*I "petscpc.h" I*/
9 
10 typedef struct {
11   Mat              fact;             /* factored matrix */
12   PetscReal        actualfill;       /* actual fill in factor */
13   PetscTruth       inplace;          /* flag indicating in-place factorization */
14   IS               row,col;          /* index sets used for reordering */
15   MatOrderingType  ordering;         /* matrix ordering */
16   PetscTruth       reuseordering;    /* reuses previous reordering computed */
17   PetscTruth       reusefill;        /* reuse fill from previous Cholesky */
18   MatFactorInfo    info;
19   MatSolverPackage solvertype;
20 } PC_Cholesky;
21 
22 EXTERN_C_BEGIN
23 #undef __FUNCT__
24 #define __FUNCT__ "PCFactorSetMatSolverPackage_Cholesky"
25 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetMatSolverPackage_Cholesky(PC pc,const MatSolverPackage stype)
26 {
27   PetscErrorCode ierr;
28   PC_Cholesky    *cholesky = (PC_Cholesky*)pc->data;
29 
30   PetscFunctionBegin;
31   if (cholesky->fact) {
32     const MatSolverPackage ltype;
33     PetscTruth             flg;
34     ierr = MatFactorGetSolverPackage(cholesky->fact,&ltype);CHKERRQ(ierr);
35     ierr = PetscStrcmp(stype,ltype,&flg);CHKERRQ(ierr);
36     if (!flg) {
37       SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Cannot change solver matrix package after PC has been setup or used");
38     }
39   }
40   ierr = PetscStrfree(cholesky->solvertype);CHKERRQ(ierr);
41   ierr = PetscStrallocpy(stype,&cholesky->solvertype);CHKERRQ(ierr);
42   PetscFunctionReturn(0);
43 }
44 EXTERN_C_END
45 
46 EXTERN_C_BEGIN
47 #undef __FUNCT__
48 #define __FUNCT__ "PCFactorSetZeroPivot_Cholesky"
49 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetZeroPivot_Cholesky(PC pc,PetscReal z)
50 {
51   PC_Cholesky *ch;
52 
53   PetscFunctionBegin;
54   ch                 = (PC_Cholesky*)pc->data;
55   ch->info.zeropivot = z;
56   PetscFunctionReturn(0);
57 }
58 EXTERN_C_END
59 
60 EXTERN_C_BEGIN
61 #undef __FUNCT__
62 #define __FUNCT__ "PCFactorSetShiftNonzero_Cholesky"
63 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftNonzero_Cholesky(PC pc,PetscReal shift)
64 {
65   PC_Cholesky *dir;
66 
67   PetscFunctionBegin;
68   dir = (PC_Cholesky*)pc->data;
69   if (shift == (PetscReal) PETSC_DECIDE) {
70     dir->info.shiftnz = 1.e-12;
71   } else {
72     dir->info.shiftnz = shift;
73   }
74   PetscFunctionReturn(0);
75 }
76 EXTERN_C_END
77 
78 EXTERN_C_BEGIN
79 #undef __FUNCT__
80 #define __FUNCT__ "PCFactorSetShiftPd_Cholesky"
81 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftPd_Cholesky(PC pc,PetscTruth shift)
82 {
83   PC_Cholesky *dir;
84 
85   PetscFunctionBegin;
86   dir = (PC_Cholesky*)pc->data;
87   if (shift) {
88     dir->info.shiftpd = 1.0;
89   } else {
90     dir->info.shiftpd = 0.0;
91   }
92   PetscFunctionReturn(0);
93 }
94 EXTERN_C_END
95 
96 EXTERN_C_BEGIN
97 #undef __FUNCT__
98 #define __FUNCT__ "PCFactorSetReuseOrdering_Cholesky"
99 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseOrdering_Cholesky(PC pc,PetscTruth flag)
100 {
101   PC_Cholesky *lu;
102 
103   PetscFunctionBegin;
104   lu               = (PC_Cholesky*)pc->data;
105   lu->reuseordering = flag;
106   PetscFunctionReturn(0);
107 }
108 EXTERN_C_END
109 
110 EXTERN_C_BEGIN
111 #undef __FUNCT__
112 #define __FUNCT__ "PCFactorSetReuseFill_Cholesky"
113 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseFill_Cholesky(PC pc,PetscTruth flag)
114 {
115   PC_Cholesky *lu;
116 
117   PetscFunctionBegin;
118   lu = (PC_Cholesky*)pc->data;
119   lu->reusefill = flag;
120   PetscFunctionReturn(0);
121 }
122 EXTERN_C_END
123 
124 #undef __FUNCT__
125 #define __FUNCT__ "PCSetFromOptions_Cholesky"
126 static PetscErrorCode PCSetFromOptions_Cholesky(PC pc)
127 {
128   PC_Cholesky    *lu = (PC_Cholesky*)pc->data;
129   PetscErrorCode ierr;
130   PetscTruth     flg;
131   char           tname[256], solvertype[64];
132   PetscFList     ordlist;
133 
134   PetscFunctionBegin;
135   ierr = MatOrderingRegisterAll(PETSC_NULL);CHKERRQ(ierr);
136   ierr = PetscOptionsHead("Cholesky options");CHKERRQ(ierr);
137     ierr = PetscOptionsName("-pc_factor_in_place","Form Cholesky in the same memory as the matrix","PCFactorSetUseInPlace",&flg);CHKERRQ(ierr);
138     if (flg) {
139       ierr = PCFactorSetUseInPlace(pc);CHKERRQ(ierr);
140     }
141     ierr = PetscOptionsReal("-pc_factor_fill","Expected non-zeros in Cholesky/non-zeros in matrix","PCFactorSetFill",lu->info.fill,&lu->info.fill,0);CHKERRQ(ierr);
142 
143     ierr = PetscOptionsName("-pc_factor_reuse_fill","Use fill from previous factorization","PCFactorSetReuseFill",&flg);CHKERRQ(ierr);
144     if (flg) {
145       ierr = PCFactorSetReuseFill(pc,PETSC_TRUE);CHKERRQ(ierr);
146     }
147     ierr = PetscOptionsName("-pc_factor_reuse_ordering","Reuse ordering from previous factorization","PCFactorSetReuseOrdering",&flg);CHKERRQ(ierr);
148     if (flg) {
149       ierr = PCFactorSetReuseOrdering(pc,PETSC_TRUE);CHKERRQ(ierr);
150     }
151 
152     ierr = MatGetOrderingList(&ordlist);CHKERRQ(ierr);
153     ierr = PetscOptionsList("-pc_factor_mat_ordering_type","Reordering to reduce nonzeros in Cholesky","PCFactorSetMatOrderingType",ordlist,lu->ordering,tname,256,&flg);CHKERRQ(ierr);
154     if (flg) {
155       ierr = PCFactorSetMatOrderingType(pc,tname);CHKERRQ(ierr);
156     }
157 
158     /* maybe should have MatGetSolverTypes(Mat,&list) like the ordering list */
159     ierr = PetscOptionsString("-pc_factor_mat_solver_package","Specific Cholesky solver to use","MatGetFactor",lu->solvertype,solvertype,64,&flg);CHKERRQ(ierr);
160     if (flg) {
161       ierr = PCFactorSetMatSolverPackage(pc,solvertype);CHKERRQ(ierr);
162     }
163 
164     ierr = PetscOptionsName("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",&flg);CHKERRQ(ierr);
165     if (flg) {
166       ierr = PCFactorSetShiftNonzero(pc,(PetscReal) PETSC_DECIDE);CHKERRQ(ierr);
167     }
168     ierr = PetscOptionsReal("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",lu->info.shiftnz,&lu->info.shiftnz,0);CHKERRQ(ierr);
169     ierr = PetscOptionsName("-pc_factor_shift_positive_definite","Manteuffel shift applied to diagonal","PCFactorSetShiftPd",&flg);CHKERRQ(ierr);
170     if (flg) {
171       ierr = PCFactorSetShiftPd(pc,PETSC_TRUE);CHKERRQ(ierr);
172     }
173     ierr = PetscOptionsReal("-pc_factor_zeropivot","Pivot is considered zero if less than","PCFactorSetZeroPivot",lu->info.zeropivot,&lu->info.zeropivot,0);CHKERRQ(ierr);
174 
175   ierr = PetscOptionsTail();CHKERRQ(ierr);
176   PetscFunctionReturn(0);
177 }
178 
179 #undef __FUNCT__
180 #define __FUNCT__ "PCView_Cholesky"
181 static PetscErrorCode PCView_Cholesky(PC pc,PetscViewer viewer)
182 {
183   PC_Cholesky    *lu = (PC_Cholesky*)pc->data;
184   PetscErrorCode ierr;
185   PetscTruth     iascii,isstring;
186 
187   PetscFunctionBegin;
188   ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr);
189   ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);CHKERRQ(ierr);
190   if (iascii) {
191 
192     if (lu->inplace) {ierr = PetscViewerASCIIPrintf(viewer,"  Cholesky: in-place factorization\n");CHKERRQ(ierr);}
193     else             {ierr = PetscViewerASCIIPrintf(viewer,"  Cholesky: out-of-place factorization\n");CHKERRQ(ierr);}
194     ierr = PetscViewerASCIIPrintf(viewer,"    matrix ordering: %s\n",lu->ordering);CHKERRQ(ierr);
195     if (lu->fact) {
196       ierr = PetscViewerASCIIPrintf(viewer,"  Cholesky: factor fill ratio needed %G\n",lu->actualfill);CHKERRQ(ierr);
197       ierr = PetscViewerASCIIPrintf(viewer,"       Factored matrix follows\n");CHKERRQ(ierr);
198       ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
199       ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
200       ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
201       ierr = PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);CHKERRQ(ierr);
202       ierr = MatView(lu->fact,viewer);CHKERRQ(ierr);
203       ierr = PetscViewerPopFormat(viewer);CHKERRQ(ierr);
204       ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
205       ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
206       ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
207     }
208     if (lu->reusefill)    {ierr = PetscViewerASCIIPrintf(viewer,"       Reusing fill from past factorization\n");CHKERRQ(ierr);}
209     if (lu->reuseordering) {ierr = PetscViewerASCIIPrintf(viewer,"       Reusing reordering from past factorization\n");CHKERRQ(ierr);}
210   } else if (isstring) {
211     ierr = PetscViewerStringSPrintf(viewer," order=%s",lu->ordering);CHKERRQ(ierr);CHKERRQ(ierr);
212   } else {
213     SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCCHOLESKY",((PetscObject)viewer)->type_name);
214   }
215   PetscFunctionReturn(0);
216 }
217 
218 #undef __FUNCT__
219 #define __FUNCT__ "PCFactorGetMatrix_Cholesky"
220 static PetscErrorCode PCFactorGetMatrix_Cholesky(PC pc,Mat *mat)
221 {
222   PC_Cholesky *dir = (PC_Cholesky*)pc->data;
223 
224   PetscFunctionBegin;
225   if (!dir->fact) SETERRQ(PETSC_ERR_ORDER,"Matrix not yet factored; call after KSPSetUp() or PCSetUp()");
226   *mat = dir->fact;
227   PetscFunctionReturn(0);
228 }
229 
230 #undef __FUNCT__
231 #define __FUNCT__ "PCSetUp_Cholesky"
232 static PetscErrorCode PCSetUp_Cholesky(PC pc)
233 {
234   PetscErrorCode ierr;
235   PetscTruth     flg;
236   PC_Cholesky    *dir = (PC_Cholesky*)pc->data;
237 
238   PetscFunctionBegin;
239   if (dir->reusefill && pc->setupcalled) dir->info.fill = dir->actualfill;
240 
241   if (dir->inplace) {
242     if (dir->row && dir->col && (dir->row != dir->col)) {
243       ierr = ISDestroy(dir->row);CHKERRQ(ierr);
244       dir->row = 0;
245     }
246     if (dir->col) {
247       ierr = ISDestroy(dir->col);CHKERRQ(ierr);
248       dir->col = 0;
249     }
250     ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr);
251     if (dir->col && (dir->row != dir->col)) {  /* only use row ordering for SBAIJ */
252       ierr = ISDestroy(dir->col);CHKERRQ(ierr);
253       dir->col=0;
254     }
255     if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);}
256     ierr = MatCholeskyFactor(pc->pmat,dir->row,&dir->info);CHKERRQ(ierr);
257     dir->fact = pc->pmat;
258   } else {
259     MatInfo info;
260     if (!pc->setupcalled) {
261       ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr);
262       /* check if dir->row == dir->col */
263       ierr = ISEqual(dir->row,dir->col,&flg);CHKERRQ(ierr);
264       if (!flg) SETERRQ(PETSC_ERR_ARG_INCOMP,"row and column permutations must equal");
265       ierr = ISDestroy(dir->col);CHKERRQ(ierr); /* only pass one ordering into CholeskyFactor */
266       dir->col=0;
267 
268       ierr = PetscOptionsHasName(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&flg);CHKERRQ(ierr);
269       if (flg) {
270         PetscReal tol = 1.e-10;
271         ierr = PetscOptionsGetReal(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&tol,PETSC_NULL);CHKERRQ(ierr);
272         ierr = MatReorderForNonzeroDiagonal(pc->pmat,tol,dir->row,dir->row);CHKERRQ(ierr);
273       }
274       if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);}
275       ierr = MatGetFactor(pc->pmat,dir->solvertype,MAT_FACTOR_CHOLESKY,&dir->fact);CHKERRQ(ierr);
276       ierr = MatCholeskyFactorSymbolic(dir->fact,pc->pmat,dir->row,&dir->info);CHKERRQ(ierr);
277       ierr = MatGetInfo(dir->fact,MAT_LOCAL,&info);CHKERRQ(ierr);
278       dir->actualfill = info.fill_ratio_needed;
279       ierr = PetscLogObjectParent(pc,dir->fact);CHKERRQ(ierr);
280     } else if (pc->flag != SAME_NONZERO_PATTERN) {
281       if (!dir->reuseordering) {
282         if (dir->row && dir->col && (dir->row != dir->col)) {
283           ierr = ISDestroy(dir->row);CHKERRQ(ierr);
284           dir->row = 0;
285         }
286         if (dir->col) {
287           ierr = ISDestroy(dir->col);CHKERRQ(ierr);
288           dir->col =0;
289         }
290         ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr);
291         if (dir->col && (dir->row != dir->col)) {  /* only use row ordering for SBAIJ */
292           ierr = ISDestroy(dir->col);CHKERRQ(ierr);
293           dir->col=0;
294         }
295         ierr = PetscOptionsHasName(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&flg);CHKERRQ(ierr);
296         if (flg) {
297           PetscReal tol = 1.e-10;
298           ierr = PetscOptionsGetReal(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&tol,PETSC_NULL);CHKERRQ(ierr);
299           ierr = MatReorderForNonzeroDiagonal(pc->pmat,tol,dir->row,dir->row);CHKERRQ(ierr);
300         }
301         if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);}
302       }
303       ierr = MatDestroy(dir->fact);CHKERRQ(ierr);
304       ierr = MatGetFactor(pc->pmat,dir->solvertype,MAT_FACTOR_CHOLESKY,&dir->fact);CHKERRQ(ierr);
305       ierr = MatCholeskyFactorSymbolic(dir->fact,pc->pmat,dir->row,&dir->info);CHKERRQ(ierr);
306       ierr = MatGetInfo(dir->fact,MAT_LOCAL,&info);CHKERRQ(ierr);
307       dir->actualfill = info.fill_ratio_needed;
308       ierr = PetscLogObjectParent(pc,dir->fact);CHKERRQ(ierr);
309     }
310     ierr = MatCholeskyFactorNumeric(dir->fact,pc->pmat,&dir->info);CHKERRQ(ierr);
311   }
312   PetscFunctionReturn(0);
313 }
314 
315 #undef __FUNCT__
316 #define __FUNCT__ "PCDestroy_Cholesky"
317 static PetscErrorCode PCDestroy_Cholesky(PC pc)
318 {
319   PC_Cholesky    *dir = (PC_Cholesky*)pc->data;
320   PetscErrorCode ierr;
321 
322   PetscFunctionBegin;
323   if (!dir->inplace && dir->fact) {ierr = MatDestroy(dir->fact);CHKERRQ(ierr);}
324   if (dir->row) {ierr = ISDestroy(dir->row);CHKERRQ(ierr);}
325   if (dir->col) {ierr = ISDestroy(dir->col);CHKERRQ(ierr);}
326   ierr = PetscStrfree(dir->ordering);CHKERRQ(ierr);
327   ierr = PetscStrfree(dir->solvertype);CHKERRQ(ierr);
328   ierr = PetscFree(dir);CHKERRQ(ierr);
329   PetscFunctionReturn(0);
330 }
331 
332 #undef __FUNCT__
333 #define __FUNCT__ "PCApply_Cholesky"
334 static PetscErrorCode PCApply_Cholesky(PC pc,Vec x,Vec y)
335 {
336   PC_Cholesky    *dir = (PC_Cholesky*)pc->data;
337   PetscErrorCode ierr;
338 
339   PetscFunctionBegin;
340   if (dir->inplace) {ierr = MatSolve(pc->pmat,x,y);CHKERRQ(ierr);}
341   else              {ierr = MatSolve(dir->fact,x,y);CHKERRQ(ierr);}
342   PetscFunctionReturn(0);
343 }
344 
345 #undef __FUNCT__
346 #define __FUNCT__ "PCApplyTranspose_Cholesky"
347 static PetscErrorCode PCApplyTranspose_Cholesky(PC pc,Vec x,Vec y)
348 {
349   PC_Cholesky    *dir = (PC_Cholesky*)pc->data;
350   PetscErrorCode ierr;
351 
352   PetscFunctionBegin;
353   if (dir->inplace) {ierr = MatSolveTranspose(pc->pmat,x,y);CHKERRQ(ierr);}
354   else              {ierr = MatSolveTranspose(dir->fact,x,y);CHKERRQ(ierr);}
355   PetscFunctionReturn(0);
356 }
357 
358 /* -----------------------------------------------------------------------------------*/
359 
360 EXTERN_C_BEGIN
361 #undef __FUNCT__
362 #define __FUNCT__ "PCFactorSetFill_Cholesky"
363 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetFill_Cholesky(PC pc,PetscReal fill)
364 {
365   PC_Cholesky *dir;
366 
367   PetscFunctionBegin;
368   dir = (PC_Cholesky*)pc->data;
369   dir->info.fill = fill;
370   PetscFunctionReturn(0);
371 }
372 EXTERN_C_END
373 
374 EXTERN_C_BEGIN
375 #undef __FUNCT__
376 #define __FUNCT__ "PCFactorSetUseInPlace_Cholesky"
377 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetUseInPlace_Cholesky(PC pc)
378 {
379   PC_Cholesky *dir;
380 
381   PetscFunctionBegin;
382   dir = (PC_Cholesky*)pc->data;
383   dir->inplace = PETSC_TRUE;
384   PetscFunctionReturn(0);
385 }
386 EXTERN_C_END
387 
388 EXTERN_C_BEGIN
389 #undef __FUNCT__
390 #define __FUNCT__ "PCFactorSetMatOrderingType_Cholesky"
391 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetMatOrderingType_Cholesky(PC pc,const MatOrderingType ordering)
392 {
393   PC_Cholesky    *dir = (PC_Cholesky*)pc->data;
394   PetscErrorCode ierr;
395 
396   PetscFunctionBegin;
397   ierr = PetscStrfree(dir->ordering);CHKERRQ(ierr);
398   ierr = PetscStrallocpy(ordering,&dir->ordering);CHKERRQ(ierr);
399   PetscFunctionReturn(0);
400 }
401 EXTERN_C_END
402 
403 /* -----------------------------------------------------------------------------------*/
404 
405 #undef __FUNCT__
406 #define __FUNCT__ "PCFactorSetReuseOrdering"
407 /*@
408    PCFactorSetReuseOrdering - When similar matrices are factored, this
409    causes the ordering computed in the first factor to be used for all
410    following factors.
411 
412    Collective on PC
413 
414    Input Parameters:
415 +  pc - the preconditioner context
416 -  flag - PETSC_TRUE to reuse else PETSC_FALSE
417 
418    Options Database Key:
419 .  -pc_factor_reuse_ordering - Activate PCFactorSetReuseOrdering()
420 
421    Level: intermediate
422 
423 .keywords: PC, levels, reordering, factorization, incomplete, LU
424 
425 .seealso: PCFactorSetReuseFill()
426 @*/
427 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseOrdering(PC pc,PetscTruth flag)
428 {
429   PetscErrorCode ierr,(*f)(PC,PetscTruth);
430 
431   PetscFunctionBegin;
432   PetscValidHeaderSpecific(pc,PC_COOKIE,1);
433   ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetReuseOrdering_C",(void (**)(void))&f);CHKERRQ(ierr);
434   if (f) {
435     ierr = (*f)(pc,flag);CHKERRQ(ierr);
436   }
437   PetscFunctionReturn(0);
438 }
439 
440 #undef __FUNCT__
441 #define __FUNCT__ "PCFactorSetReuseFill"
442 /*@
443    PCFactorSetReuseFill - When matrices with same different nonzero structure are factored,
444    this causes later ones to use the fill ratio computed in the initial factorization.
445 
446    Collective on PC
447 
448    Input Parameters:
449 +  pc - the preconditioner context
450 -  flag - PETSC_TRUE to reuse else PETSC_FALSE
451 
452    Options Database Key:
453 .  -pc_factor_reuse_fill - Activates PCFactorSetReuseFill()
454 
455    Level: intermediate
456 
457 .keywords: PC, levels, reordering, factorization, incomplete, Cholesky
458 
459 .seealso: PCFactorSetReuseOrdering()
460 @*/
461 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseFill(PC pc,PetscTruth flag)
462 {
463   PetscErrorCode ierr,(*f)(PC,PetscTruth);
464 
465   PetscFunctionBegin;
466   PetscValidHeaderSpecific(pc,PC_COOKIE,2);
467   ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetReuseFill_C",(void (**)(void))&f);CHKERRQ(ierr);
468   if (f) {
469     ierr = (*f)(pc,flag);CHKERRQ(ierr);
470   }
471   PetscFunctionReturn(0);
472 }
473 
474 /*MC
475    PCCHOLESKY - Uses a direct solver, based on Cholesky factorization, as a preconditioner
476 
477    Options Database Keys:
478 +  -pc_factor_reuse_ordering - Activate PCFactorSetReuseOrdering()
479 .  -pc_factor_mat_solver_package - Actives PCFactorSetMatSolverPackage() to choose the direct solver, like spooles
480 .  -pc_factor_reuse_fill - Activates PCFactorSetReuseFill()
481 .  -pc_factor_fill <fill> - Sets fill amount
482 .  -pc_factor_in_place - Activates in-place factorization
483 .  -pc_factor_mat_ordering_type <nd,rcm,...> - Sets ordering routine
484 .  -pc_factor_shift_nonzero <shift> - Sets shift amount or PETSC_DECIDE for the default
485 -  -pc_factor_shift_positive_definite [PETSC_TRUE/PETSC_FALSE] - Activate/Deactivate PCFactorSetShiftPd(); the value
486    is optional with PETSC_TRUE being the default
487 
488    Notes: Not all options work for all matrix formats
489 
490    Level: beginner
491 
492    Concepts: Cholesky factorization, direct solver
493 
494    Notes: Usually this will compute an "exact" solution in one iteration and does
495           not need a Krylov method (i.e. you can use -ksp_type preonly, or
496           KSPSetType(ksp,KSPPREONLY) for the Krylov method
497 
498 .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PC,
499            PCILU, PCLU, PCICC, PCFactorSetReuseOrdering(), PCFactorSetReuseFill(), PCFactorGetMatrix(),
500            PCFactorSetFill(), PCFactorSetShiftNonzero(), PCFactorSetShiftPd(),
501 	   PCFactorSetUseInPlace(), PCFactorSetMatOrderingType(),PCFactorSetShiftNonzero(),PCFactorSetShiftPd()
502 
503 M*/
504 
505 EXTERN_C_BEGIN
506 #undef __FUNCT__
507 #define __FUNCT__ "PCCreate_Cholesky"
508 PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_Cholesky(PC pc)
509 {
510   PetscErrorCode ierr;
511   PC_Cholesky    *dir;
512 
513   PetscFunctionBegin;
514   ierr = PetscNewLog(pc,PC_Cholesky,&dir);CHKERRQ(ierr);
515 
516   dir->fact                   = 0;
517   dir->inplace                = PETSC_FALSE;
518   ierr = MatFactorInfoInitialize(&dir->info);CHKERRQ(ierr);
519   dir->info.fill              = 5.0;
520   dir->info.shiftnz           = 0.0;
521   dir->info.shiftpd           = 0.0; /* false */
522   dir->info.pivotinblocks     = 1.0;
523   dir->col                    = 0;
524   dir->row                    = 0;
525   ierr = PetscStrallocpy(MATORDERING_NATURAL,&dir->ordering);CHKERRQ(ierr);
526   ierr = PetscStrallocpy(MAT_SOLVER_PETSC,&dir->solvertype);CHKERRQ(ierr);
527   dir->reusefill        = PETSC_FALSE;
528   dir->reuseordering    = PETSC_FALSE;
529   pc->data              = (void*)dir;
530 
531   pc->ops->destroy           = PCDestroy_Cholesky;
532   pc->ops->apply             = PCApply_Cholesky;
533   pc->ops->applytranspose    = PCApplyTranspose_Cholesky;
534   pc->ops->setup             = PCSetUp_Cholesky;
535   pc->ops->setfromoptions    = PCSetFromOptions_Cholesky;
536   pc->ops->view              = PCView_Cholesky;
537   pc->ops->applyrichardson   = 0;
538   pc->ops->getfactoredmatrix = PCFactorGetMatrix_Cholesky;
539 
540   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatSolverPackage_C","PCFactorSetMatSolverPackage_Cholesky",
541                     PCFactorSetMatSolverPackage_Cholesky);CHKERRQ(ierr);
542   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetZeroPivot_C","PCFactorSetZeroPivot_Cholesky",
543                     PCFactorSetZeroPivot_Cholesky);CHKERRQ(ierr);
544   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftNonzero_C","PCFactorSetShiftNonzero_Cholesky",
545                     PCFactorSetShiftNonzero_Cholesky);CHKERRQ(ierr);
546   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftPd_C","PCFactorSetShiftPd_Cholesky",
547                     PCFactorSetShiftPd_Cholesky);CHKERRQ(ierr);
548 
549   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetFill_C","PCFactorSetFill_Cholesky",
550                     PCFactorSetFill_Cholesky);CHKERRQ(ierr);
551   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetUseInPlace_C","PCFactorSetUseInPlace_Cholesky",
552                     PCFactorSetUseInPlace_Cholesky);CHKERRQ(ierr);
553   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatOrderingType_C","PCFactorSetMatOrderingType_Cholesky",
554                     PCFactorSetMatOrderingType_Cholesky);CHKERRQ(ierr);
555   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseOrdering_C","PCFactorSetReuseOrdering_Cholesky",
556                     PCFactorSetReuseOrdering_Cholesky);CHKERRQ(ierr);
557   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseFill_C","PCFactorSetReuseFill_Cholesky",
558                     PCFactorSetReuseFill_Cholesky);CHKERRQ(ierr);
559   PetscFunctionReturn(0);
560 }
561 EXTERN_C_END
562