Actual source code: factor.c
petsc-3.4.5 2014-06-29
2: #include <../src/ksp/pc/impls/factor/factor.h> /*I "petscpc.h" I*/
6: /*@
7: PCFactorSetUpMatSolverPackage - Can be called after KSPSetOperators() or PCSetOperators(), causes MatGetFactor() to be called so then one may
8: set the options for that particular factorization object.
10: Input Parameter:
11: . pc - the preconditioner context
13: 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.
15: .seealso: PCFactorSetMatSolverPackage(), PCFactorGetMatrix()
17: Level: intermediate
19: @*/
20: PetscErrorCode PCFactorSetUpMatSolverPackage(PC pc)
21: {
26: PetscTryMethod(pc,"PCFactorSetUpMatSolverPackage_C",(PC),(pc));
27: return(0);
28: }
32: /*@
33: PCFactorSetZeroPivot - Sets the size at which smaller pivots are declared to be zero
35: Logically Collective on PC
37: Input Parameters:
38: + pc - the preconditioner context
39: - zero - all pivots smaller than this will be considered zero
41: Options Database Key:
42: . -pc_factor_zeropivot <zero> - Sets tolerance for what is considered a zero pivot
44: Level: intermediate
46: .keywords: PC, set, factorization, direct, fill
48: .seealso: PCFactorSetShiftType(), PCFactorSetShiftAmount()
49: @*/
50: PetscErrorCode PCFactorSetZeroPivot(PC pc,PetscReal zero)
51: {
57: PetscTryMethod(pc,"PCFactorSetZeroPivot_C",(PC,PetscReal),(pc,zero));
58: return(0);
59: }
63: /*@
64: PCFactorSetShiftType - adds a particular type of quantity to the diagonal of the matrix during
65: numerical factorization, thus the matrix has nonzero pivots
67: Logically Collective on PC
69: Input Parameters:
70: + pc - the preconditioner context
71: - shifttype - type of shift; one of MAT_SHIFT_NONE, MAT_SHIFT_NONZERO, MAT_SHIFT_POSITIVE_DEFINITE, MAT_SHIFT_INBLOCKS
73: Options Database Key:
74: . -pc_factor_shift_type <shifttype> - Sets shift type or PETSC_DECIDE for the default; use '-help' for a list of available types
76: Level: intermediate
78: .keywords: PC, set, factorization,
80: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftAmount()
81: @*/
82: PetscErrorCode PCFactorSetShiftType(PC pc,MatFactorShiftType shifttype)
83: {
89: PetscTryMethod(pc,"PCFactorSetShiftType_C",(PC,MatFactorShiftType),(pc,shifttype));
90: return(0);
91: }
95: /*@
96: PCFactorSetShiftAmount - adds a quantity to the diagonal of the matrix during
97: numerical factorization, thus the matrix has nonzero pivots
99: Logically Collective on PC
101: Input Parameters:
102: + pc - the preconditioner context
103: - shiftamount - amount of shift
105: Options Database Key:
106: . -pc_factor_shift_amount <shiftamount> - Sets shift amount or PETSC_DECIDE for the default
108: Level: intermediate
110: .keywords: PC, set, factorization,
112: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftType()
113: @*/
114: PetscErrorCode PCFactorSetShiftAmount(PC pc,PetscReal shiftamount)
115: {
121: PetscTryMethod(pc,"PCFactorSetShiftAmount_C",(PC,PetscReal),(pc,shiftamount));
122: return(0);
123: }
127: /*
128: PCFactorSetDropTolerance - The preconditioner will use an ILU
129: based on a drop tolerance. (Under development)
131: Logically Collective on PC
133: Input Parameters:
134: + pc - the preconditioner context
135: . dt - the drop tolerance, try from 1.e-10 to .1
136: . dtcol - tolerance for column pivot, good values [0.1 to 0.01]
137: - maxrowcount - the max number of nonzeros allowed in a row, best value
138: depends on the number of nonzeros in row of original matrix
140: Options Database Key:
141: . -pc_factor_drop_tolerance <dt,dtcol,maxrowcount> - Sets drop tolerance
143: Level: intermediate
145: There are NO default values for the 3 parameters, you must set them with reasonable values for your
146: matrix. We don't know how to compute reasonable values.
148: .keywords: PC, levels, reordering, factorization, incomplete, ILU
149: */
150: PetscErrorCode PCFactorSetDropTolerance(PC pc,PetscReal dt,PetscReal dtcol,PetscInt maxrowcount)
151: {
158: PetscTryMethod(pc,"PCFactorSetDropTolerance_C",(PC,PetscReal,PetscReal,PetscInt),(pc,dt,dtcol,maxrowcount));
159: return(0);
160: }
164: /*@
165: PCFactorGetLevels - Gets the number of levels of fill to use.
167: Logically Collective on PC
169: Input Parameters:
170: . pc - the preconditioner context
172: Output Parameter:
173: . levels - number of levels of fill
175: Level: intermediate
177: .keywords: PC, levels, fill, factorization, incomplete, ILU
178: @*/
179: PetscErrorCode PCFactorGetLevels(PC pc,PetscInt *levels)
180: {
185: PetscTryMethod(pc,"PCFactorGetLevels_C",(PC,PetscInt*),(pc,levels));
186: return(0);
187: }
191: /*@
192: PCFactorSetLevels - Sets the number of levels of fill to use.
194: Logically Collective on PC
196: Input Parameters:
197: + pc - the preconditioner context
198: - levels - number of levels of fill
200: Options Database Key:
201: . -pc_factor_levels <levels> - Sets fill level
203: Level: intermediate
205: .keywords: PC, levels, fill, factorization, incomplete, ILU
206: @*/
207: PetscErrorCode PCFactorSetLevels(PC pc,PetscInt levels)
208: {
213: if (levels < 0) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"negative levels");
215: PetscTryMethod(pc,"PCFactorSetLevels_C",(PC,PetscInt),(pc,levels));
216: return(0);
217: }
221: /*@
222: PCFactorSetAllowDiagonalFill - Causes all diagonal matrix entries to be
223: treated as level 0 fill even if there is no non-zero location.
225: Logically Collective on PC
227: Input Parameters:
228: + pc - the preconditioner context
230: Options Database Key:
231: . -pc_factor_diagonal_fill
233: Notes:
234: Does not apply with 0 fill.
236: Level: intermediate
238: .keywords: PC, levels, fill, factorization, incomplete, ILU
239: @*/
240: PetscErrorCode PCFactorSetAllowDiagonalFill(PC pc)
241: {
246: PetscTryMethod(pc,"PCFactorSetAllowDiagonalFill_C",(PC),(pc));
247: return(0);
248: }
252: /*@
253: PCFactorReorderForNonzeroDiagonal - reorders rows/columns of matrix to remove zeros from diagonal
255: Logically Collective on PC
257: Input Parameters:
258: + pc - the preconditioner context
259: - tol - diagonal entries smaller than this in absolute value are considered zero
261: Options Database Key:
262: . -pc_factor_nonzeros_along_diagonal
264: Level: intermediate
266: .keywords: PC, set, factorization, direct, fill
268: .seealso: PCFactorSetFill(), PCFactorSetShiftNonzero(), PCFactorSetZeroPivot(), MatReorderForNonzeroDiagonal()
269: @*/
270: PetscErrorCode PCFactorReorderForNonzeroDiagonal(PC pc,PetscReal rtol)
271: {
277: PetscTryMethod(pc,"PCFactorReorderForNonzeroDiagonal_C",(PC,PetscReal),(pc,rtol));
278: return(0);
279: }
283: /*@C
284: PCFactorSetMatSolverPackage - sets the software that is used to perform the factorization
286: Logically Collective on PC
288: Input Parameters:
289: + pc - the preconditioner context
290: - stype - for example, superlu, superlu_dist
292: Options Database Key:
293: . -pc_factor_mat_solver_package <stype> - petsc, superlu, superlu_dist, mumps, cusparse
295: Level: intermediate
297: Note:
298: By default this will use the PETSc factorization if it exists
301: .keywords: PC, set, factorization, direct, fill
303: .seealso: MatGetFactor(), MatSolverPackage, PCFactorGetMatSolverPackage()
305: @*/
306: PetscErrorCode PCFactorSetMatSolverPackage(PC pc,const MatSolverPackage stype)
307: {
312: PetscTryMethod(pc,"PCFactorSetMatSolverPackage_C",(PC,const MatSolverPackage),(pc,stype));
313: return(0);
314: }
318: /*@C
319: PCFactorGetMatSolverPackage - gets the software that is used to perform the factorization
321: Not Collective
323: Input Parameter:
324: . pc - the preconditioner context
326: Output Parameter:
327: . stype - for example, superlu, superlu_dist (NULL if the PC does not have a solver package)
329: Level: intermediate
332: .keywords: PC, set, factorization, direct, fill
334: .seealso: MatGetFactor(), MatSolverPackage, PCFactorGetMatSolverPackage()
336: @*/
337: PetscErrorCode PCFactorGetMatSolverPackage(PC pc,const MatSolverPackage *stype)
338: {
339: PetscErrorCode ierr,(*f)(PC,const MatSolverPackage*);
343: PetscObjectQueryFunction((PetscObject)pc,"PCFactorGetMatSolverPackage_C",&f);
344: if (f) {
345: (*f)(pc,stype);
346: } else {
347: *stype = NULL;
348: }
349: return(0);
350: }
354: /*@
355: PCFactorSetFill - Indicate the amount of fill you expect in the factored matrix,
356: fill = number nonzeros in factor/number nonzeros in original matrix.
358: Not Collective, each process can expect a different amount of fill
360: Input Parameters:
361: + pc - the preconditioner context
362: - fill - amount of expected fill
364: Options Database Key:
365: . -pc_factor_fill <fill> - Sets fill amount
367: Level: intermediate
369: Note:
370: For sparse matrix factorizations it is difficult to predict how much
371: fill to expect. By running with the option -info PETSc will print the
372: actual amount of fill used; allowing you to set the value accurately for
373: future runs. Default PETSc uses a value of 5.0
375: This parameter has NOTHING to do with the levels-of-fill of ILU(). That is set with PCFactorSetLevels() or -pc_factor_levels.
378: .keywords: PC, set, factorization, direct, fill
380: @*/
381: PetscErrorCode PCFactorSetFill(PC pc,PetscReal fill)
382: {
387: if (fill < 1.0) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Fill factor cannot be less then 1.0");
388: PetscTryMethod(pc,"PCFactorSetFill_C",(PC,PetscReal),(pc,fill));
389: return(0);
390: }
394: /*@
395: PCFactorSetUseInPlace - Tells the system to do an in-place factorization.
396: For dense matrices, this enables the solution of much larger problems.
397: For sparse matrices the factorization cannot be done truly in-place
398: so this does not save memory during the factorization, but after the matrix
399: is factored, the original unfactored matrix is freed, thus recovering that
400: space.
402: Logically Collective on PC
404: Input Parameters:
405: . pc - the preconditioner context
407: Options Database Key:
408: . -pc_factor_in_place - Activates in-place factorization
410: Notes:
411: PCFactorSetUseInplace() can only be used with the KSP method KSPPREONLY or when
412: a different matrix is provided for the multiply and the preconditioner in
413: a call to KSPSetOperators().
414: This is because the Krylov space methods require an application of the
415: matrix multiplication, which is not possible here because the matrix has
416: been factored in-place, replacing the original matrix.
418: Level: intermediate
420: .keywords: PC, set, factorization, direct, inplace, in-place, LU
422: .seealso: PCILUSetUseInPlace()
423: @*/
424: PetscErrorCode PCFactorSetUseInPlace(PC pc)
425: {
430: PetscTryMethod(pc,"PCFactorSetUseInPlace_C",(PC),(pc));
431: return(0);
432: }
436: /*@C
437: PCFactorSetMatOrderingType - Sets the ordering routine (to reduce fill) to
438: be used in the LU factorization.
440: Logically Collective on PC
442: Input Parameters:
443: + pc - the preconditioner context
444: - ordering - the matrix ordering name, for example, MATORDERINGND or MATORDERINGRCM
446: Options Database Key:
447: . -pc_factor_mat_ordering_type <nd,rcm,...> - Sets ordering routine
449: Level: intermediate
451: Notes: nested dissection is used by default
453: For Cholesky and ICC and the SBAIJ format reorderings are not available,
454: since only the upper triangular part of the matrix is stored. You can use the
455: SeqAIJ format in this case to get reorderings.
457: @*/
458: PetscErrorCode PCFactorSetMatOrderingType(PC pc,MatOrderingType ordering)
459: {
464: PetscTryMethod(pc,"PCFactorSetMatOrderingType_C",(PC,MatOrderingType),(pc,ordering));
465: return(0);
466: }
470: /*@
471: PCFactorSetColumnPivot - Determines when column pivoting is done during matrix factorization.
472: For PETSc dense matrices column pivoting is always done, for PETSc sparse matrices
473: it is never done. For the MATLAB and SuperLU factorization this is used.
475: Logically Collective on PC
477: Input Parameters:
478: + pc - the preconditioner context
479: - dtcol - 0.0 implies no pivoting, 1.0 complete pivoting (slower, requires more memory but more stable)
481: Options Database Key:
482: . -pc_factor_pivoting <dtcol>
484: Level: intermediate
486: .seealso: PCILUSetMatOrdering(), PCFactorSetPivotInBlocks()
487: @*/
488: PetscErrorCode PCFactorSetColumnPivot(PC pc,PetscReal dtcol)
489: {
495: PetscTryMethod(pc,"PCFactorSetColumnPivot_C",(PC,PetscReal),(pc,dtcol));
496: return(0);
497: }
501: /*@
502: PCFactorSetPivotInBlocks - Determines if pivoting is done while factoring each block
503: with BAIJ or SBAIJ matrices
505: Logically Collective on PC
507: Input Parameters:
508: + pc - the preconditioner context
509: - pivot - PETSC_TRUE or PETSC_FALSE
511: Options Database Key:
512: . -pc_factor_pivot_in_blocks <true,false>
514: Level: intermediate
516: .seealso: PCILUSetMatOrdering(), PCFactorSetColumnPivot()
517: @*/
518: PetscErrorCode PCFactorSetPivotInBlocks(PC pc,PetscBool pivot)
519: {
525: PetscTryMethod(pc,"PCFactorSetPivotInBlocks_C",(PC,PetscBool),(pc,pivot));
526: return(0);
527: }
531: /*@
532: PCFactorSetReuseFill - When matrices with same different nonzero structure are factored,
533: this causes later ones to use the fill ratio computed in the initial factorization.
535: Logically Collective on PC
537: Input Parameters:
538: + pc - the preconditioner context
539: - flag - PETSC_TRUE to reuse else PETSC_FALSE
541: Options Database Key:
542: . -pc_factor_reuse_fill - Activates PCFactorSetReuseFill()
544: Level: intermediate
546: .keywords: PC, levels, reordering, factorization, incomplete, Cholesky
548: .seealso: PCFactorSetReuseOrdering()
549: @*/
550: PetscErrorCode PCFactorSetReuseFill(PC pc,PetscBool flag)
551: {
557: PetscTryMethod(pc,"PCFactorSetReuseFill_C",(PC,PetscBool),(pc,flag));
558: return(0);
559: }