Actual source code: factor.c
petsc-3.8.4 2018-03-24
2: #include <../src/ksp/pc/impls/factor/factor.h>
4: static PetscErrorCode PCFactorSetReuseOrdering_Factor(PC pc,PetscBool flag)
5: {
6: PC_Factor *lu = (PC_Factor*)pc->data;
9: lu->reuseordering = flag;
10: return(0);
11: }
13: static PetscErrorCode PCFactorSetReuseFill_Factor(PC pc,PetscBool flag)
14: {
15: PC_Factor *lu = (PC_Factor*)pc->data;
18: lu->reusefill = flag;
19: return(0);
20: }
22: static PetscErrorCode PCFactorSetUseInPlace_Factor(PC pc,PetscBool flg)
23: {
24: PC_Factor *dir = (PC_Factor*)pc->data;
27: dir->inplace = flg;
28: return(0);
29: }
31: static PetscErrorCode PCFactorGetUseInPlace_Factor(PC pc,PetscBool *flg)
32: {
33: PC_Factor *dir = (PC_Factor*)pc->data;
36: *flg = dir->inplace;
37: return(0);
38: }
40: /*@
41: PCFactorSetUpMatSolverPackage - Can be called after KSPSetOperators() or PCSetOperators(), causes MatGetFactor() to be called so then one may
42: set the options for that particular factorization object.
44: Input Parameter:
45: . pc - the preconditioner context
47: Notes: After you have called this function (which has to be after the KSPSetOperators() or PCSetOperators()) you can call PCFactorGetMatrix() and then set factor options on that matrix.
49: .seealso: PCFactorSetMatSolverPackage(), PCFactorGetMatrix()
51: Level: intermediate
53: @*/
54: PetscErrorCode PCFactorSetUpMatSolverPackage(PC pc)
55: {
60: PetscTryMethod(pc,"PCFactorSetUpMatSolverPackage_C",(PC),(pc));
61: return(0);
62: }
64: /*@
65: PCFactorSetZeroPivot - Sets the size at which smaller pivots are declared to be zero
67: Logically Collective on PC
69: Input Parameters:
70: + pc - the preconditioner context
71: - zero - all pivots smaller than this will be considered zero
73: Options Database Key:
74: . -pc_factor_zeropivot <zero> - Sets tolerance for what is considered a zero pivot
76: Level: intermediate
78: .keywords: PC, set, factorization, direct, fill
80: .seealso: PCFactorSetShiftType(), PCFactorSetShiftAmount()
81: @*/
82: PetscErrorCode PCFactorSetZeroPivot(PC pc,PetscReal zero)
83: {
89: PetscTryMethod(pc,"PCFactorSetZeroPivot_C",(PC,PetscReal),(pc,zero));
90: return(0);
91: }
93: /*@
94: PCFactorSetShiftType - adds a particular type of quantity to the diagonal of the matrix during
95: numerical factorization, thus the matrix has nonzero pivots
97: Logically Collective on PC
99: Input Parameters:
100: + pc - the preconditioner context
101: - shifttype - type of shift; one of MAT_SHIFT_NONE, MAT_SHIFT_NONZERO, MAT_SHIFT_POSITIVE_DEFINITE, MAT_SHIFT_INBLOCKS
103: Options Database Key:
104: . -pc_factor_shift_type <shifttype> - Sets shift type or PETSC_DECIDE for the default; use '-help' for a list of available types
106: Level: intermediate
108: .keywords: PC, set, factorization,
110: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftAmount()
111: @*/
112: PetscErrorCode PCFactorSetShiftType(PC pc,MatFactorShiftType shifttype)
113: {
119: PetscTryMethod(pc,"PCFactorSetShiftType_C",(PC,MatFactorShiftType),(pc,shifttype));
120: return(0);
121: }
123: /*@
124: PCFactorSetShiftAmount - adds a quantity to the diagonal of the matrix during
125: numerical factorization, thus the matrix has nonzero pivots
127: Logically Collective on PC
129: Input Parameters:
130: + pc - the preconditioner context
131: - shiftamount - amount of shift
133: Options Database Key:
134: . -pc_factor_shift_amount <shiftamount> - Sets shift amount or PETSC_DECIDE for the default
136: Level: intermediate
138: .keywords: PC, set, factorization,
140: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftType()
141: @*/
142: PetscErrorCode PCFactorSetShiftAmount(PC pc,PetscReal shiftamount)
143: {
149: PetscTryMethod(pc,"PCFactorSetShiftAmount_C",(PC,PetscReal),(pc,shiftamount));
150: return(0);
151: }
153: /*
154: PCFactorSetDropTolerance - The preconditioner will use an ILU
155: based on a drop tolerance. (Under development)
157: Logically Collective on PC
159: Input Parameters:
160: + pc - the preconditioner context
161: . dt - the drop tolerance, try from 1.e-10 to .1
162: . dtcol - tolerance for column pivot, good values [0.1 to 0.01]
163: - maxrowcount - the max number of nonzeros allowed in a row, best value
164: depends on the number of nonzeros in row of original matrix
166: Options Database Key:
167: . -pc_factor_drop_tolerance <dt,dtcol,maxrowcount> - Sets drop tolerance
169: Level: intermediate
171: There are NO default values for the 3 parameters, you must set them with reasonable values for your
172: matrix. We don't know how to compute reasonable values.
174: .keywords: PC, levels, reordering, factorization, incomplete, ILU
175: */
176: PetscErrorCode PCFactorSetDropTolerance(PC pc,PetscReal dt,PetscReal dtcol,PetscInt maxrowcount)
177: {
184: PetscTryMethod(pc,"PCFactorSetDropTolerance_C",(PC,PetscReal,PetscReal,PetscInt),(pc,dt,dtcol,maxrowcount));
185: return(0);
186: }
188: /*@
189: PCFactorGetZeroPivot - Gets the tolerance used to define a zero privot
191: Not Collective
193: Input Parameters:
194: . pc - the preconditioner context
196: Output Parameter:
197: . pivot - the tolerance
199: Level: intermediate
202: .seealso: PCFactorSetZeroPivot()
203: @*/
204: PetscErrorCode PCFactorGetZeroPivot(PC pc,PetscReal *pivot)
205: {
210: PetscUseMethod(pc,"PCFactorGetZeroPivot_C",(PC,PetscReal*),(pc,pivot));
211: return(0);
212: }
214: /*@
215: PCFactorGetShiftAmount - Gets the tolerance used to define a zero privot
217: Not Collective
219: Input Parameters:
220: . pc - the preconditioner context
222: Output Parameter:
223: . shift - how much to shift the diagonal entry
225: Level: intermediate
228: .seealso: PCFactorSetShiftAmount(), PCFactorSetShiftType(), PCFactorGetShiftType()
229: @*/
230: PetscErrorCode PCFactorGetShiftAmount(PC pc,PetscReal *shift)
231: {
236: PetscUseMethod(pc,"PCFactorGetShiftAmount_C",(PC,PetscReal*),(pc,shift));
237: return(0);
238: }
240: /*@
241: PCFactorGetShiftType - Gets the type of shift, if any, done when a zero pivot is detected
243: Not Collective
245: Input Parameters:
246: . pc - the preconditioner context
248: Output Parameter:
249: . type - one of MAT_SHIFT_NONE, MAT_SHIFT_NONZERO, MAT_SHIFT_POSITIVE_DEFINITE, or MAT_SHIFT_INBLOCKS
251: Level: intermediate
254: .seealso: PCFactorSetShiftType(), MatFactorShiftType, PCFactorSetShiftAmount(), PCFactorGetShiftAmount()
255: @*/
256: PetscErrorCode PCFactorGetShiftType(PC pc,MatFactorShiftType *type)
257: {
262: PetscUseMethod(pc,"PCFactorGetShiftType_C",(PC,MatFactorShiftType*),(pc,type));
263: return(0);
264: }
266: /*@
267: PCFactorGetLevels - Gets the number of levels of fill to use.
269: Logically Collective on PC
271: Input Parameters:
272: . pc - the preconditioner context
274: Output Parameter:
275: . levels - number of levels of fill
277: Level: intermediate
279: .keywords: PC, levels, fill, factorization, incomplete, ILU
280: @*/
281: PetscErrorCode PCFactorGetLevels(PC pc,PetscInt *levels)
282: {
287: PetscUseMethod(pc,"PCFactorGetLevels_C",(PC,PetscInt*),(pc,levels));
288: return(0);
289: }
291: /*@
292: PCFactorSetLevels - Sets the number of levels of fill to use.
294: Logically Collective on PC
296: Input Parameters:
297: + pc - the preconditioner context
298: - levels - number of levels of fill
300: Options Database Key:
301: . -pc_factor_levels <levels> - Sets fill level
303: Level: intermediate
305: .keywords: PC, levels, fill, factorization, incomplete, ILU
306: @*/
307: PetscErrorCode PCFactorSetLevels(PC pc,PetscInt levels)
308: {
313: if (levels < 0) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"negative levels");
315: PetscTryMethod(pc,"PCFactorSetLevels_C",(PC,PetscInt),(pc,levels));
316: return(0);
317: }
319: /*@
320: PCFactorSetAllowDiagonalFill - Causes all diagonal matrix entries to be
321: treated as level 0 fill even if there is no non-zero location.
323: Logically Collective on PC
325: Input Parameters:
326: + pc - the preconditioner context
327: - flg - PETSC_TRUE to turn on, PETSC_FALSE to turn off
329: Options Database Key:
330: . -pc_factor_diagonal_fill
332: Notes:
333: Does not apply with 0 fill.
335: Level: intermediate
337: .keywords: PC, levels, fill, factorization, incomplete, ILU
339: .seealso: PCFactorGetAllowDiagonalFill()
341: @*/
342: PetscErrorCode PCFactorSetAllowDiagonalFill(PC pc,PetscBool flg)
343: {
348: PetscTryMethod(pc,"PCFactorSetAllowDiagonalFill_C",(PC,PetscBool),(pc,flg));
349: return(0);
350: }
352: /*@
353: PCFactorGetAllowDiagonalFill - Determines if all diagonal matrix entries are
354: treated as level 0 fill even if there is no non-zero location.
356: Logically Collective on PC
358: Input Parameter:
359: . pc - the preconditioner context
361: Output Parameter:
362: . flg - PETSC_TRUE to turn on, PETSC_FALSE to turn off
364: Options Database Key:
365: . -pc_factor_diagonal_fill
367: Notes:
368: Does not apply with 0 fill.
370: Level: intermediate
372: .keywords: PC, levels, fill, factorization, incomplete, ILU
374: .seealso: PCFactorSetAllowDiagonalFill()
376: @*/
377: PetscErrorCode PCFactorGetAllowDiagonalFill(PC pc,PetscBool *flg)
378: {
383: PetscUseMethod(pc,"PCFactorGetAllowDiagonalFill_C",(PC,PetscBool*),(pc,flg));
384: return(0);
385: }
387: /*@
388: PCFactorReorderForNonzeroDiagonal - reorders rows/columns of matrix to remove zeros from diagonal
390: Logically Collective on PC
392: Input Parameters:
393: + pc - the preconditioner context
394: - tol - diagonal entries smaller than this in absolute value are considered zero
396: Options Database Key:
397: . -pc_factor_nonzeros_along_diagonal <tol>
399: Level: intermediate
401: .keywords: PC, set, factorization, direct, fill
403: .seealso: PCFactorSetFill(), PCFactorSetShiftNonzero(), PCFactorSetZeroPivot(), MatReorderForNonzeroDiagonal()
404: @*/
405: PetscErrorCode PCFactorReorderForNonzeroDiagonal(PC pc,PetscReal rtol)
406: {
412: PetscTryMethod(pc,"PCFactorReorderForNonzeroDiagonal_C",(PC,PetscReal),(pc,rtol));
413: return(0);
414: }
416: /*@C
417: PCFactorSetMatSolverPackage - sets the software that is used to perform the factorization
419: Logically Collective on PC
421: Input Parameters:
422: + pc - the preconditioner context
423: - stype - for example, superlu, superlu_dist
425: Options Database Key:
426: . -pc_factor_mat_solver_package <stype> - petsc, superlu, superlu_dist, mumps, cusparse
428: Level: intermediate
430: Note:
431: By default this will use the PETSc factorization if it exists
434: .keywords: PC, set, factorization, direct, fill
436: .seealso: MatGetFactor(), MatSolverPackage, PCFactorGetMatSolverPackage()
438: @*/
439: PetscErrorCode PCFactorSetMatSolverPackage(PC pc,const MatSolverPackage stype)
440: {
445: PetscTryMethod(pc,"PCFactorSetMatSolverPackage_C",(PC,const MatSolverPackage),(pc,stype));
446: return(0);
447: }
449: /*@C
450: PCFactorGetMatSolverPackage - gets the software that is used to perform the factorization
452: Not Collective
454: Input Parameter:
455: . pc - the preconditioner context
457: Output Parameter:
458: . stype - for example, superlu, superlu_dist (NULL if the PC does not have a solver package)
460: Level: intermediate
463: .keywords: PC, set, factorization, direct, fill
465: .seealso: MatGetFactor(), MatSolverPackage, PCFactorGetMatSolverPackage()
467: @*/
468: PetscErrorCode PCFactorGetMatSolverPackage(PC pc,const MatSolverPackage *stype)
469: {
470: PetscErrorCode ierr,(*f)(PC,const MatSolverPackage*);
474: PetscObjectQueryFunction((PetscObject)pc,"PCFactorGetMatSolverPackage_C",&f);
475: if (f) {
476: (*f)(pc,stype);
477: } else {
478: *stype = NULL;
479: }
480: return(0);
481: }
483: /*@
484: PCFactorSetFill - Indicate the amount of fill you expect in the factored matrix,
485: fill = number nonzeros in factor/number nonzeros in original matrix.
487: Not Collective, each process can expect a different amount of fill
489: Input Parameters:
490: + pc - the preconditioner context
491: - fill - amount of expected fill
493: Options Database Key:
494: . -pc_factor_fill <fill> - Sets fill amount
496: Level: intermediate
498: Note:
499: For sparse matrix factorizations it is difficult to predict how much
500: fill to expect. By running with the option -info PETSc will print the
501: actual amount of fill used; allowing you to set the value accurately for
502: future runs. Default PETSc uses a value of 5.0
504: This parameter has NOTHING to do with the levels-of-fill of ILU(). That is set with PCFactorSetLevels() or -pc_factor_levels.
507: .keywords: PC, set, factorization, direct, fill
509: @*/
510: PetscErrorCode PCFactorSetFill(PC pc,PetscReal fill)
511: {
516: if (fill < 1.0) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Fill factor cannot be less then 1.0");
517: PetscTryMethod(pc,"PCFactorSetFill_C",(PC,PetscReal),(pc,fill));
518: return(0);
519: }
521: /*@
522: PCFactorSetUseInPlace - Tells the system to do an in-place factorization.
523: For dense matrices, this enables the solution of much larger problems.
524: For sparse matrices the factorization cannot be done truly in-place
525: so this does not save memory during the factorization, but after the matrix
526: is factored, the original unfactored matrix is freed, thus recovering that
527: space.
529: Logically Collective on PC
531: Input Parameters:
532: + pc - the preconditioner context
533: - flg - PETSC_TRUE to enable, PETSC_FALSE to disable
535: Options Database Key:
536: . -pc_factor_in_place <true,false>- Activate/deactivate in-place factorization
538: Notes:
539: PCFactorSetUseInplace() can only be used with the KSP method KSPPREONLY or when
540: a different matrix is provided for the multiply and the preconditioner in
541: a call to KSPSetOperators().
542: This is because the Krylov space methods require an application of the
543: matrix multiplication, which is not possible here because the matrix has
544: been factored in-place, replacing the original matrix.
546: Level: intermediate
548: .keywords: PC, set, factorization, direct, inplace, in-place, LU
550: .seealso: PCFactorGetUseInPlace()
551: @*/
552: PetscErrorCode PCFactorSetUseInPlace(PC pc,PetscBool flg)
553: {
558: PetscTryMethod(pc,"PCFactorSetUseInPlace_C",(PC,PetscBool),(pc,flg));
559: return(0);
560: }
562: /*@
563: PCFactorGetUseInPlace - Determines if an in-place factorization is being used.
565: Logically Collective on PC
567: Input Parameter:
568: . pc - the preconditioner context
570: Output Parameter:
571: . flg - PETSC_TRUE to enable, PETSC_FALSE to disable
573: Level: intermediate
575: .keywords: PC, set, factorization, direct, inplace, in-place, LU
577: .seealso: PCFactorSetUseInPlace()
578: @*/
579: PetscErrorCode PCFactorGetUseInPlace(PC pc,PetscBool *flg)
580: {
585: PetscUseMethod(pc,"PCFactorGetUseInPlace_C",(PC,PetscBool*),(pc,flg));
586: return(0);
587: }
589: /*@C
590: PCFactorSetMatOrderingType - Sets the ordering routine (to reduce fill) to
591: be used in the LU factorization.
593: Logically Collective on PC
595: Input Parameters:
596: + pc - the preconditioner context
597: - ordering - the matrix ordering name, for example, MATORDERINGND or MATORDERINGRCM
599: Options Database Key:
600: . -pc_factor_mat_ordering_type <nd,rcm,...> - Sets ordering routine
602: Level: intermediate
604: Notes: nested dissection is used by default
606: For Cholesky and ICC and the SBAIJ format reorderings are not available,
607: since only the upper triangular part of the matrix is stored. You can use the
608: SeqAIJ format in this case to get reorderings.
610: @*/
611: PetscErrorCode PCFactorSetMatOrderingType(PC pc,MatOrderingType ordering)
612: {
617: PetscTryMethod(pc,"PCFactorSetMatOrderingType_C",(PC,MatOrderingType),(pc,ordering));
618: return(0);
619: }
621: /*@
622: PCFactorSetColumnPivot - Determines when column pivoting is done during matrix factorization.
623: For PETSc dense matrices column pivoting is always done, for PETSc sparse matrices
624: it is never done. For the MATLAB and SuperLU factorization this is used.
626: Logically Collective on PC
628: Input Parameters:
629: + pc - the preconditioner context
630: - dtcol - 0.0 implies no pivoting, 1.0 complete pivoting (slower, requires more memory but more stable)
632: Options Database Key:
633: . -pc_factor_pivoting <dtcol>
635: Level: intermediate
637: .seealso: PCILUSetMatOrdering(), PCFactorSetPivotInBlocks()
638: @*/
639: PetscErrorCode PCFactorSetColumnPivot(PC pc,PetscReal dtcol)
640: {
646: PetscTryMethod(pc,"PCFactorSetColumnPivot_C",(PC,PetscReal),(pc,dtcol));
647: return(0);
648: }
650: /*@
651: PCFactorSetPivotInBlocks - Determines if pivoting is done while factoring each block
652: with BAIJ or SBAIJ matrices
654: Logically Collective on PC
656: Input Parameters:
657: + pc - the preconditioner context
658: - pivot - PETSC_TRUE or PETSC_FALSE
660: Options Database Key:
661: . -pc_factor_pivot_in_blocks <true,false>
663: Level: intermediate
665: .seealso: PCILUSetMatOrdering(), PCFactorSetColumnPivot()
666: @*/
667: PetscErrorCode PCFactorSetPivotInBlocks(PC pc,PetscBool pivot)
668: {
674: PetscTryMethod(pc,"PCFactorSetPivotInBlocks_C",(PC,PetscBool),(pc,pivot));
675: return(0);
676: }
678: /*@
679: PCFactorSetReuseFill - When matrices with different nonzero structure are factored,
680: this causes later ones to use the fill ratio computed in the initial factorization.
682: Logically Collective on PC
684: Input Parameters:
685: + pc - the preconditioner context
686: - flag - PETSC_TRUE to reuse else PETSC_FALSE
688: Options Database Key:
689: . -pc_factor_reuse_fill - Activates PCFactorSetReuseFill()
691: Level: intermediate
693: .keywords: PC, levels, reordering, factorization, incomplete, Cholesky
695: .seealso: PCFactorSetReuseOrdering()
696: @*/
697: PetscErrorCode PCFactorSetReuseFill(PC pc,PetscBool flag)
698: {
704: PetscTryMethod(pc,"PCFactorSetReuseFill_C",(PC,PetscBool),(pc,flag));
705: return(0);
706: }
708: PetscErrorCode PCFactorInitialize(PC pc)
709: {
711: PC_Factor *fact = (PC_Factor*)pc->data;
714: MatFactorInfoInitialize(&fact->info);
715: fact->info.shifttype = (PetscReal)MAT_SHIFT_NONE;
716: fact->info.shiftamount = 100.0*PETSC_MACHINE_EPSILON;
717: fact->info.zeropivot = 100.0*PETSC_MACHINE_EPSILON;
718: fact->info.pivotinblocks = 1.0;
719: pc->ops->getfactoredmatrix = PCFactorGetMatrix_Factor;
721: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetZeroPivot_C",PCFactorSetZeroPivot_Factor);
722: PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetZeroPivot_C",PCFactorGetZeroPivot_Factor);
723: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftType_C",PCFactorSetShiftType_Factor);
724: PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetShiftType_C",PCFactorGetShiftType_Factor);
725: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftAmount_C",PCFactorSetShiftAmount_Factor);
726: PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetShiftAmount_C",PCFactorGetShiftAmount_Factor);
727: PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetMatSolverPackage_C",PCFactorGetMatSolverPackage_Factor);
728: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatSolverPackage_C",PCFactorSetMatSolverPackage_Factor);
729: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUpMatSolverPackage_C",PCFactorSetUpMatSolverPackage_Factor);
730: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetFill_C",PCFactorSetFill_Factor);
731: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatOrderingType_C",PCFactorSetMatOrderingType_Factor);
732: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetLevels_C",PCFactorSetLevels_Factor);
733: PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetLevels_C",PCFactorGetLevels_Factor);
734: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetAllowDiagonalFill_C",PCFactorSetAllowDiagonalFill_Factor);
735: PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetAllowDiagonalFill_C",PCFactorGetAllowDiagonalFill_Factor);
736: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetPivotInBlocks_C",PCFactorSetPivotInBlocks_Factor);
737: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUseInPlace_C",PCFactorSetUseInPlace_Factor);
738: PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetUseInPlace_C",PCFactorGetUseInPlace_Factor);
739: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseOrdering_C",PCFactorSetReuseOrdering_Factor);
740: PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseFill_C",PCFactorSetReuseFill_Factor);
741: return(0);
742: }