Actual source code: factimpl.c
2: #include <../src/ksp/pc/impls/factor/factor.h>
4: PetscErrorCode PCFactorSetUpMatSolverType_Factor(PC pc)
5: {
6: PC_Factor *icc = (PC_Factor *)pc->data;
9: if (!pc->setupcalled && !((PC_Factor *)icc)->fact) MatGetFactor(pc->pmat, ((PC_Factor *)icc)->solvertype, ((PC_Factor *)icc)->factortype, &((PC_Factor *)icc)->fact);
10: return 0;
11: }
13: PetscErrorCode PCFactorSetZeroPivot_Factor(PC pc, PetscReal z)
14: {
15: PC_Factor *ilu = (PC_Factor *)pc->data;
17: ilu->info.zeropivot = z;
18: return 0;
19: }
21: PetscErrorCode PCFactorSetShiftType_Factor(PC pc, MatFactorShiftType shifttype)
22: {
23: PC_Factor *dir = (PC_Factor *)pc->data;
25: if (shifttype == (MatFactorShiftType)PETSC_DECIDE) dir->info.shifttype = (PetscReal)MAT_SHIFT_NONE;
26: else {
27: dir->info.shifttype = (PetscReal)shifttype;
28: if ((shifttype == MAT_SHIFT_NONZERO || shifttype == MAT_SHIFT_INBLOCKS) && dir->info.shiftamount == 0.0) { dir->info.shiftamount = 100.0 * PETSC_MACHINE_EPSILON; /* set default amount if user has not called PCFactorSetShiftAmount() yet */ }
29: }
30: return 0;
31: }
33: PetscErrorCode PCFactorSetShiftAmount_Factor(PC pc, PetscReal shiftamount)
34: {
35: PC_Factor *dir = (PC_Factor *)pc->data;
37: if (shiftamount == (PetscReal)PETSC_DECIDE) dir->info.shiftamount = 100.0 * PETSC_MACHINE_EPSILON;
38: else dir->info.shiftamount = shiftamount;
39: return 0;
40: }
42: PetscErrorCode PCFactorSetDropTolerance_Factor(PC pc, PetscReal dt, PetscReal dtcol, PetscInt dtcount)
43: {
44: PC_Factor *ilu = (PC_Factor *)pc->data;
46: if (pc->setupcalled && (!ilu->info.usedt || ((PC_Factor *)ilu)->info.dt != dt || ((PC_Factor *)ilu)->info.dtcol != dtcol || ((PC_Factor *)ilu)->info.dtcount != dtcount)) {
47: SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "Cannot change tolerance after use");
48: }
49: ilu->info.usedt = PETSC_TRUE;
50: ilu->info.dt = dt;
51: ilu->info.dtcol = dtcol;
52: ilu->info.dtcount = dtcount;
53: ilu->info.fill = PETSC_DEFAULT;
54: return 0;
55: }
57: PetscErrorCode PCFactorSetFill_Factor(PC pc, PetscReal fill)
58: {
59: PC_Factor *dir = (PC_Factor *)pc->data;
61: dir->info.fill = fill;
62: return 0;
63: }
65: PetscErrorCode PCFactorSetMatOrderingType_Factor(PC pc, MatOrderingType ordering)
66: {
67: PC_Factor *dir = (PC_Factor *)pc->data;
68: PetscBool flg;
70: if (!pc->setupcalled) {
71: PetscFree(dir->ordering);
72: PetscStrallocpy(ordering, (char **)&dir->ordering);
73: } else {
74: PetscStrcmp(dir->ordering, ordering, &flg);
76: }
77: return 0;
78: }
80: PetscErrorCode PCFactorGetLevels_Factor(PC pc, PetscInt *levels)
81: {
82: PC_Factor *ilu = (PC_Factor *)pc->data;
84: *levels = ilu->info.levels;
85: return 0;
86: }
88: PetscErrorCode PCFactorGetZeroPivot_Factor(PC pc, PetscReal *pivot)
89: {
90: PC_Factor *ilu = (PC_Factor *)pc->data;
92: *pivot = ilu->info.zeropivot;
93: return 0;
94: }
96: PetscErrorCode PCFactorGetShiftAmount_Factor(PC pc, PetscReal *shift)
97: {
98: PC_Factor *ilu = (PC_Factor *)pc->data;
100: *shift = ilu->info.shiftamount;
101: return 0;
102: }
104: PetscErrorCode PCFactorGetShiftType_Factor(PC pc, MatFactorShiftType *type)
105: {
106: PC_Factor *ilu = (PC_Factor *)pc->data;
108: *type = (MatFactorShiftType)(int)ilu->info.shifttype;
109: return 0;
110: }
112: PetscErrorCode PCFactorSetLevels_Factor(PC pc, PetscInt levels)
113: {
114: PC_Factor *ilu = (PC_Factor *)pc->data;
116: if (!pc->setupcalled) ilu->info.levels = levels;
117: else if (ilu->info.levels != levels) {
118: PetscUseTypeMethod(pc, reset); /* remove previous factored matrices */
119: pc->setupcalled = 0; /* force a complete rebuild of preconditioner factored matrices */
120: ilu->info.levels = levels;
122: return 0;
123: }
125: PetscErrorCode PCFactorSetAllowDiagonalFill_Factor(PC pc, PetscBool flg)
126: {
127: PC_Factor *dir = (PC_Factor *)pc->data;
129: dir->info.diagonal_fill = (PetscReal)flg;
130: return 0;
131: }
133: PetscErrorCode PCFactorGetAllowDiagonalFill_Factor(PC pc, PetscBool *flg)
134: {
135: PC_Factor *dir = (PC_Factor *)pc->data;
137: *flg = dir->info.diagonal_fill ? PETSC_TRUE : PETSC_FALSE;
138: return 0;
139: }
141: PetscErrorCode PCFactorSetPivotInBlocks_Factor(PC pc, PetscBool pivot)
142: {
143: PC_Factor *dir = (PC_Factor *)pc->data;
145: dir->info.pivotinblocks = pivot ? 1.0 : 0.0;
146: return 0;
147: }
149: PetscErrorCode PCFactorGetMatrix_Factor(PC pc, Mat *mat)
150: {
151: PC_Factor *ilu = (PC_Factor *)pc->data;
154: *mat = ilu->fact;
155: return 0;
156: }
158: /* allow access to preallocation information */
159: #include <petsc/private/matimpl.h>
161: PetscErrorCode PCFactorSetMatSolverType_Factor(PC pc, MatSolverType stype)
162: {
163: PC_Factor *lu = (PC_Factor *)pc->data;
165: if (lu->fact && lu->fact->assembled) {
166: MatSolverType ltype;
167: PetscBool flg;
168: MatFactorGetSolverType(lu->fact, <ype);
169: PetscStrcmp(stype, ltype, &flg);
171: }
173: PetscFree(lu->solvertype);
174: PetscStrallocpy(stype, &lu->solvertype);
175: return 0;
176: }
178: PetscErrorCode PCFactorGetMatSolverType_Factor(PC pc, MatSolverType *stype)
179: {
180: PC_Factor *lu = (PC_Factor *)pc->data;
182: *stype = lu->solvertype;
183: return 0;
184: }
186: PetscErrorCode PCFactorSetColumnPivot_Factor(PC pc, PetscReal dtcol)
187: {
188: PC_Factor *dir = (PC_Factor *)pc->data;
191: dir->info.dtcol = dtcol;
192: return 0;
193: }
195: PetscErrorCode PCSetFromOptions_Factor(PC pc, PetscOptionItems *PetscOptionsObject)
196: {
197: PC_Factor *factor = (PC_Factor *)pc->data;
198: PetscBool flg, set;
199: char tname[256], solvertype[64];
200: PetscFunctionList ordlist;
201: PetscEnum etmp;
202: PetscBool inplace;
204: PCFactorGetUseInPlace(pc, &inplace);
205: PetscOptionsBool("-pc_factor_in_place", "Form factored matrix in the same memory as the matrix", "PCFactorSetUseInPlace", inplace, &flg, &set);
206: if (set) PCFactorSetUseInPlace(pc, flg);
207: PetscOptionsReal("-pc_factor_fill", "Expected non-zeros in factored matrix", "PCFactorSetFill", ((PC_Factor *)factor)->info.fill, &((PC_Factor *)factor)->info.fill, NULL);
209: PetscOptionsEnum("-pc_factor_shift_type", "Type of shift to add to diagonal", "PCFactorSetShiftType", MatFactorShiftTypes, (PetscEnum)(int)((PC_Factor *)factor)->info.shifttype, &etmp, &flg);
210: if (flg) PCFactorSetShiftType(pc, (MatFactorShiftType)etmp);
211: PetscOptionsReal("-pc_factor_shift_amount", "Shift added to diagonal", "PCFactorSetShiftAmount", ((PC_Factor *)factor)->info.shiftamount, &((PC_Factor *)factor)->info.shiftamount, NULL);
213: PetscOptionsReal("-pc_factor_zeropivot", "Pivot is considered zero if less than", "PCFactorSetZeroPivot", ((PC_Factor *)factor)->info.zeropivot, &((PC_Factor *)factor)->info.zeropivot, NULL);
214: 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);
216: 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);
217: if (set) PCFactorSetPivotInBlocks(pc, flg);
219: PetscOptionsBool("-pc_factor_reuse_fill", "Use fill from previous factorization", "PCFactorSetReuseFill", PETSC_FALSE, &flg, &set);
220: if (set) PCFactorSetReuseFill(pc, flg);
221: PetscOptionsBool("-pc_factor_reuse_ordering", "Reuse ordering from previous factorization", "PCFactorSetReuseOrdering", PETSC_FALSE, &flg, &set);
222: if (set) PCFactorSetReuseOrdering(pc, flg);
224: PetscOptionsDeprecated("-pc_factor_mat_solver_package", "-pc_factor_mat_solver_type", "3.9", NULL);
225: PetscOptionsString("-pc_factor_mat_solver_type", "Specific direct solver to use", "MatGetFactor", ((PC_Factor *)factor)->solvertype, solvertype, sizeof(solvertype), &flg);
226: if (flg) PCFactorSetMatSolverType(pc, solvertype);
227: PCFactorSetDefaultOrdering_Factor(pc);
228: MatGetOrderingList(&ordlist);
229: PetscOptionsFList("-pc_factor_mat_ordering_type", "Reordering to reduce nonzeros in factored matrix", "PCFactorSetMatOrderingType", ordlist, ((PC_Factor *)factor)->ordering, tname, sizeof(tname), &flg);
230: if (flg) PCFactorSetMatOrderingType(pc, tname);
231: return 0;
232: }
234: PetscErrorCode PCView_Factor(PC pc, PetscViewer viewer)
235: {
236: PC_Factor *factor = (PC_Factor *)pc->data;
237: PetscBool isstring, iascii, canuseordering;
238: MatInfo info;
239: MatOrderingType ordering;
240: PetscViewerFormat format;
242: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSTRING, &isstring);
243: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
244: if (iascii) {
245: if (factor->inplace) {
246: PetscViewerASCIIPrintf(viewer, " in-place factorization\n");
247: } else {
248: PetscViewerASCIIPrintf(viewer, " out-of-place factorization\n");
249: }
251: if (factor->reusefill) PetscViewerASCIIPrintf(viewer, " Reusing fill from past factorization\n");
252: if (factor->reuseordering) PetscViewerASCIIPrintf(viewer, " Reusing reordering from past factorization\n");
253: if (factor->factortype == MAT_FACTOR_ILU || factor->factortype == MAT_FACTOR_ICC) {
254: if (factor->info.dt > 0) {
255: PetscViewerASCIIPrintf(viewer, " drop tolerance %g\n", (double)factor->info.dt);
256: PetscViewerASCIIPrintf(viewer, " max nonzeros per row %" PetscInt_FMT "\n", (PetscInt)factor->info.dtcount);
257: PetscViewerASCIIPrintf(viewer, " column permutation tolerance %g\n", (double)factor->info.dtcol);
258: } else if (factor->info.levels == 1) {
259: PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " level of fill\n", (PetscInt)factor->info.levels);
260: } else {
261: PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " levels of fill\n", (PetscInt)factor->info.levels);
262: }
263: }
265: PetscViewerASCIIPrintf(viewer, " tolerance for zero pivot %g\n", (double)factor->info.zeropivot);
266: if (MatFactorShiftTypesDetail[(int)factor->info.shifttype]) { /* Only print when using a nontrivial shift */
267: PetscViewerASCIIPrintf(viewer, " using %s [%s]\n", MatFactorShiftTypesDetail[(int)factor->info.shifttype], MatFactorShiftTypes[(int)factor->info.shifttype]);
268: }
270: if (factor->fact) {
271: MatFactorGetCanUseOrdering(factor->fact, &canuseordering);
272: if (!canuseordering) ordering = MATORDERINGEXTERNAL;
273: else ordering = factor->ordering;
274: PetscViewerASCIIPrintf(viewer, " matrix ordering: %s\n", ordering);
275: if (!factor->fact->assembled) {
276: PetscViewerASCIIPrintf(viewer, " matrix solver type: %s\n", factor->fact->solvertype);
277: PetscViewerASCIIPrintf(viewer, " matrix not yet factored; no additional information available\n");
278: } else {
279: MatGetInfo(factor->fact, MAT_LOCAL, &info);
280: PetscViewerASCIIPrintf(viewer, " factor fill ratio given %g, needed %g\n", (double)info.fill_ratio_given, (double)info.fill_ratio_needed);
281: PetscViewerASCIIPrintf(viewer, " Factored matrix follows:\n");
282: PetscViewerASCIIPushTab(viewer);
283: PetscViewerASCIIPushTab(viewer);
284: PetscViewerASCIIPushTab(viewer);
285: PetscViewerGetFormat(viewer, &format);
286: PetscViewerPushFormat(viewer, format != PETSC_VIEWER_ASCII_INFO_DETAIL ? PETSC_VIEWER_ASCII_INFO : PETSC_VIEWER_ASCII_INFO_DETAIL);
287: MatView(factor->fact, viewer);
288: PetscViewerPopFormat(viewer);
289: PetscViewerASCIIPopTab(viewer);
290: PetscViewerASCIIPopTab(viewer);
291: PetscViewerASCIIPopTab(viewer);
292: }
293: }
295: } else if (isstring) {
296: MatFactorType t;
297: MatGetFactorType(factor->fact, &t);
298: if (t == MAT_FACTOR_ILU || t == MAT_FACTOR_ICC) PetscViewerStringSPrintf(viewer, " lvls=%" PetscInt_FMT ",order=%s", (PetscInt)factor->info.levels, factor->ordering);
299: }
300: return 0;
301: }