Actual source code: factimpl.c
2: #include <../src/ksp/pc/impls/factor/factor.h>
4: /* ------------------------------------------------------------------------------------------*/
6: PetscErrorCode PCFactorSetUpMatSolverType_Factor(PC pc)
7: {
8: PC_Factor *icc = (PC_Factor*)pc->data;
11: if (!pc->setupcalled && !((PC_Factor*)icc)->fact) {
12: MatGetFactor(pc->pmat,((PC_Factor*)icc)->solvertype,((PC_Factor*)icc)->factortype,&((PC_Factor*)icc)->fact);
13: }
14: return 0;
15: }
17: PetscErrorCode PCFactorSetZeroPivot_Factor(PC pc,PetscReal z)
18: {
19: PC_Factor *ilu = (PC_Factor*)pc->data;
21: ilu->info.zeropivot = z;
22: return 0;
23: }
25: PetscErrorCode PCFactorSetShiftType_Factor(PC pc,MatFactorShiftType shifttype)
26: {
27: PC_Factor *dir = (PC_Factor*)pc->data;
29: if (shifttype == (MatFactorShiftType)PETSC_DECIDE) dir->info.shifttype = (PetscReal) MAT_SHIFT_NONE;
30: else {
31: dir->info.shifttype = (PetscReal) shifttype;
32: if ((shifttype == MAT_SHIFT_NONZERO || shifttype == MAT_SHIFT_INBLOCKS) && dir->info.shiftamount == 0.0) {
33: dir->info.shiftamount = 100.0*PETSC_MACHINE_EPSILON; /* set default amount if user has not called PCFactorSetShiftAmount() yet */
34: }
35: }
36: return 0;
37: }
39: PetscErrorCode PCFactorSetShiftAmount_Factor(PC pc,PetscReal shiftamount)
40: {
41: PC_Factor *dir = (PC_Factor*)pc->data;
43: if (shiftamount == (PetscReal) PETSC_DECIDE) dir->info.shiftamount = 100.0*PETSC_MACHINE_EPSILON;
44: else dir->info.shiftamount = shiftamount;
45: return 0;
46: }
48: PetscErrorCode PCFactorSetDropTolerance_Factor(PC pc,PetscReal dt,PetscReal dtcol,PetscInt dtcount)
49: {
50: PC_Factor *ilu = (PC_Factor*)pc->data;
52: if (pc->setupcalled && (!ilu->info.usedt || ((PC_Factor*)ilu)->info.dt != dt || ((PC_Factor*)ilu)->info.dtcol != dtcol || ((PC_Factor*)ilu)->info.dtcount != dtcount)) {
53: SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Cannot change tolerance after use");
54: }
55: ilu->info.usedt = PETSC_TRUE;
56: ilu->info.dt = dt;
57: ilu->info.dtcol = dtcol;
58: ilu->info.dtcount = dtcount;
59: ilu->info.fill = PETSC_DEFAULT;
60: return 0;
61: }
63: PetscErrorCode PCFactorSetFill_Factor(PC pc,PetscReal fill)
64: {
65: PC_Factor *dir = (PC_Factor*)pc->data;
67: dir->info.fill = fill;
68: return 0;
69: }
71: PetscErrorCode PCFactorSetMatOrderingType_Factor(PC pc,MatOrderingType ordering)
72: {
73: PC_Factor *dir = (PC_Factor*)pc->data;
74: PetscBool flg;
76: if (!pc->setupcalled) {
77: PetscFree(dir->ordering);
78: PetscStrallocpy(ordering,(char**)&dir->ordering);
79: } else {
80: PetscStrcmp(dir->ordering,ordering,&flg);
82: }
83: return 0;
84: }
86: PetscErrorCode PCFactorGetLevels_Factor(PC pc,PetscInt *levels)
87: {
88: PC_Factor *ilu = (PC_Factor*)pc->data;
90: *levels = ilu->info.levels;
91: return 0;
92: }
94: PetscErrorCode PCFactorGetZeroPivot_Factor(PC pc,PetscReal *pivot)
95: {
96: PC_Factor *ilu = (PC_Factor*)pc->data;
98: *pivot = ilu->info.zeropivot;
99: return 0;
100: }
102: PetscErrorCode PCFactorGetShiftAmount_Factor(PC pc,PetscReal *shift)
103: {
104: PC_Factor *ilu = (PC_Factor*)pc->data;
106: *shift = ilu->info.shiftamount;
107: return 0;
108: }
110: PetscErrorCode PCFactorGetShiftType_Factor(PC pc,MatFactorShiftType *type)
111: {
112: PC_Factor *ilu = (PC_Factor*)pc->data;
114: *type = (MatFactorShiftType) (int) ilu->info.shifttype;
115: return 0;
116: }
118: PetscErrorCode PCFactorSetLevels_Factor(PC pc,PetscInt levels)
119: {
120: PC_Factor *ilu = (PC_Factor*)pc->data;
122: if (!pc->setupcalled) ilu->info.levels = levels;
123: else if (ilu->info.levels != levels) {
124: (*pc->ops->reset)(pc); /* remove previous factored matrices */
125: pc->setupcalled = 0; /* force a complete rebuild of preconditioner factored matrices */
126: ilu->info.levels = levels;
128: return 0;
129: }
131: PetscErrorCode PCFactorSetAllowDiagonalFill_Factor(PC pc,PetscBool flg)
132: {
133: PC_Factor *dir = (PC_Factor*)pc->data;
135: dir->info.diagonal_fill = (PetscReal) flg;
136: return 0;
137: }
139: PetscErrorCode PCFactorGetAllowDiagonalFill_Factor(PC pc,PetscBool *flg)
140: {
141: PC_Factor *dir = (PC_Factor*)pc->data;
143: *flg = dir->info.diagonal_fill ? PETSC_TRUE : PETSC_FALSE;
144: return 0;
145: }
147: /* ------------------------------------------------------------------------------------------*/
149: PetscErrorCode PCFactorSetPivotInBlocks_Factor(PC pc,PetscBool pivot)
150: {
151: PC_Factor *dir = (PC_Factor*)pc->data;
153: dir->info.pivotinblocks = pivot ? 1.0 : 0.0;
154: return 0;
155: }
157: PetscErrorCode PCFactorGetMatrix_Factor(PC pc,Mat *mat)
158: {
159: PC_Factor *ilu = (PC_Factor*)pc->data;
162: *mat = ilu->fact;
163: return 0;
164: }
166: /* allow access to preallocation information */
167: #include <petsc/private/matimpl.h>
169: PetscErrorCode PCFactorSetMatSolverType_Factor(PC pc,MatSolverType stype)
170: {
171: PC_Factor *lu = (PC_Factor*)pc->data;
173: if (lu->fact && lu->fact->assembled) {
174: MatSolverType ltype;
175: PetscBool flg;
176: MatFactorGetSolverType(lu->fact,<ype);
177: PetscStrcmp(stype,ltype,&flg);
179: }
181: PetscFree(lu->solvertype);
182: PetscStrallocpy(stype,&lu->solvertype);
183: return 0;
184: }
186: PetscErrorCode PCFactorGetMatSolverType_Factor(PC pc,MatSolverType *stype)
187: {
188: PC_Factor *lu = (PC_Factor*)pc->data;
190: *stype = lu->solvertype;
191: return 0;
192: }
194: PetscErrorCode PCFactorSetColumnPivot_Factor(PC pc,PetscReal dtcol)
195: {
196: PC_Factor *dir = (PC_Factor*)pc->data;
199: dir->info.dtcol = dtcol;
200: return 0;
201: }
203: PetscErrorCode PCSetFromOptions_Factor(PetscOptionItems *PetscOptionsObject,PC pc)
204: {
205: PC_Factor *factor = (PC_Factor*)pc->data;
206: PetscBool flg,set;
207: char tname[256], solvertype[64];
208: PetscFunctionList ordlist;
209: PetscEnum etmp;
210: PetscBool inplace;
212: PCFactorGetUseInPlace(pc,&inplace);
213: PetscOptionsBool("-pc_factor_in_place","Form factored matrix in the same memory as the matrix","PCFactorSetUseInPlace",inplace,&flg,&set);
214: if (set) {
215: PCFactorSetUseInPlace(pc,flg);
216: }
217: PetscOptionsReal("-pc_factor_fill","Expected non-zeros in factored matrix","PCFactorSetFill",((PC_Factor*)factor)->info.fill,&((PC_Factor*)factor)->info.fill,NULL);
219: PetscOptionsEnum("-pc_factor_shift_type","Type of shift to add to diagonal","PCFactorSetShiftType",MatFactorShiftTypes,(PetscEnum)(int)((PC_Factor*)factor)->info.shifttype,&etmp,&flg);
220: if (flg) {
221: PCFactorSetShiftType(pc,(MatFactorShiftType)etmp);
222: }
223: PetscOptionsReal("-pc_factor_shift_amount","Shift added to diagonal","PCFactorSetShiftAmount",((PC_Factor*)factor)->info.shiftamount,&((PC_Factor*)factor)->info.shiftamount,NULL);
225: PetscOptionsReal("-pc_factor_zeropivot","Pivot is considered zero if less than","PCFactorSetZeroPivot",((PC_Factor*)factor)->info.zeropivot,&((PC_Factor*)factor)->info.zeropivot,NULL);
226: PetscOptionsReal("-pc_factor_column_pivot","Column pivot tolerance (used only for some factorization)","PCFactorSetColumnPivot",((PC_Factor*)factor)->info.dtcol,&((PC_Factor*)factor)->info.dtcol,&flg);
228: PetscOptionsBool("-pc_factor_pivot_in_blocks","Pivot inside matrix dense blocks for BAIJ and SBAIJ","PCFactorSetPivotInBlocks",((PC_Factor*)factor)->info.pivotinblocks ? PETSC_TRUE : PETSC_FALSE,&flg,&set);
229: if (set) {
230: PCFactorSetPivotInBlocks(pc,flg);
231: }
233: PetscOptionsBool("-pc_factor_reuse_fill","Use fill from previous factorization","PCFactorSetReuseFill",PETSC_FALSE,&flg,&set);
234: if (set) {
235: PCFactorSetReuseFill(pc,flg);
236: }
237: PetscOptionsBool("-pc_factor_reuse_ordering","Reuse ordering from previous factorization","PCFactorSetReuseOrdering",PETSC_FALSE,&flg,&set);
238: if (set) {
239: PCFactorSetReuseOrdering(pc,flg);
240: }
242: PetscOptionsDeprecated("-pc_factor_mat_solver_package","-pc_factor_mat_solver_type","3.9",NULL);
243: PetscOptionsString("-pc_factor_mat_solver_type","Specific direct solver to use","MatGetFactor",((PC_Factor*)factor)->solvertype,solvertype,sizeof(solvertype),&flg);
244: if (flg) {
245: PCFactorSetMatSolverType(pc,solvertype);
246: }
247: PCFactorSetDefaultOrdering_Factor(pc);
248: MatGetOrderingList(&ordlist);
249: PetscOptionsFList("-pc_factor_mat_ordering_type","Reordering to reduce nonzeros in factored matrix","PCFactorSetMatOrderingType",ordlist,((PC_Factor*)factor)->ordering,tname,sizeof(tname),&flg);
250: if (flg) {
251: PCFactorSetMatOrderingType(pc,tname);
252: }
253: return 0;
254: }
256: PetscErrorCode PCView_Factor(PC pc,PetscViewer viewer)
257: {
258: PC_Factor *factor = (PC_Factor*)pc->data;
259: PetscBool isstring,iascii,canuseordering;
260: MatInfo info;
261: MatOrderingType ordering;
263: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);
264: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
265: if (iascii) {
266: if (factor->inplace) {
267: PetscViewerASCIIPrintf(viewer," in-place factorization\n");
268: } else {
269: PetscViewerASCIIPrintf(viewer," out-of-place factorization\n");
270: }
272: if (factor->reusefill) PetscViewerASCIIPrintf(viewer," Reusing fill from past factorization\n");
273: if (factor->reuseordering) PetscViewerASCIIPrintf(viewer," Reusing reordering from past factorization\n");
274: if (factor->factortype == MAT_FACTOR_ILU || factor->factortype == MAT_FACTOR_ICC) {
275: if (factor->info.dt > 0) {
276: PetscViewerASCIIPrintf(viewer," drop tolerance %g\n",(double)factor->info.dt);
277: PetscViewerASCIIPrintf(viewer," max nonzeros per row %D\n",factor->info.dtcount);
278: PetscViewerASCIIPrintf(viewer," column permutation tolerance %g\n",(double)factor->info.dtcol);
279: } else if (factor->info.levels == 1) {
280: PetscViewerASCIIPrintf(viewer," %D level of fill\n",(PetscInt)factor->info.levels);
281: } else {
282: PetscViewerASCIIPrintf(viewer," %D levels of fill\n",(PetscInt)factor->info.levels);
283: }
284: }
286: PetscViewerASCIIPrintf(viewer," tolerance for zero pivot %g\n",(double)factor->info.zeropivot);
287: if (MatFactorShiftTypesDetail[(int)factor->info.shifttype]) { /* Only print when using a nontrivial shift */
288: PetscViewerASCIIPrintf(viewer," using %s [%s]\n",MatFactorShiftTypesDetail[(int)factor->info.shifttype],MatFactorShiftTypes[(int)factor->info.shifttype]);
289: }
291: if (factor->fact) {
292: MatFactorGetCanUseOrdering(factor->fact,&canuseordering);
293: if (!canuseordering) ordering = MATORDERINGEXTERNAL;
294: else ordering = factor->ordering;
295: PetscViewerASCIIPrintf(viewer," matrix ordering: %s\n",ordering);
296: if (!factor->fact->assembled) {
297: PetscViewerASCIIPrintf(viewer," matrix solver type: %s\n",factor->fact->solvertype);
298: PetscViewerASCIIPrintf(viewer," matrix not yet factored; no additional information available\n");
299: } else {
300: MatGetInfo(factor->fact,MAT_LOCAL,&info);
301: PetscViewerASCIIPrintf(viewer," factor fill ratio given %g, needed %g\n",(double)info.fill_ratio_given,(double)info.fill_ratio_needed);
302: PetscViewerASCIIPrintf(viewer," Factored matrix follows:\n");
303: PetscViewerASCIIPushTab(viewer);
304: PetscViewerASCIIPushTab(viewer);
305: PetscViewerASCIIPushTab(viewer);
306: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
307: MatView(factor->fact,viewer);
308: PetscViewerPopFormat(viewer);
309: PetscViewerASCIIPopTab(viewer);
310: PetscViewerASCIIPopTab(viewer);
311: PetscViewerASCIIPopTab(viewer);
312: }
313: }
315: } else if (isstring) {
316: MatFactorType t;
317: MatGetFactorType(factor->fact,&t);
318: if (t == MAT_FACTOR_ILU || t == MAT_FACTOR_ICC) {
319: PetscViewerStringSPrintf(viewer," lvls=%D,order=%s",(PetscInt)factor->info.levels,factor->ordering);
320: }
321: }
322: return 0;
323: }