Actual source code: factor.c

petsc-master 2019-06-15
Report Typos and Errors

  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:     PCFactorSetUpMatSolverType - 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:
 48:     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.

 50: .seealso: PCFactorSetMatSolverType(), PCFactorGetMatrix()

 52:   Level: intermediate

 54: @*/
 55: PetscErrorCode PCFactorSetUpMatSolverType(PC pc)
 56: {

 61:   PetscTryMethod(pc,"PCFactorSetUpMatSolverType_C",(PC),(pc));
 62:   return(0);
 63: }

 65: /*@
 66:    PCFactorSetZeroPivot - Sets the size at which smaller pivots are declared to be zero

 68:    Logically Collective on PC

 70:    Input Parameters:
 71: +  pc - the preconditioner context
 72: -  zero - all pivots smaller than this will be considered zero

 74:    Options Database Key:
 75: .  -pc_factor_zeropivot <zero> - Sets tolerance for what is considered a zero pivot

 77:    Level: intermediate

 79: .seealso: PCFactorSetShiftType(), PCFactorSetShiftAmount()
 80: @*/
 81: PetscErrorCode  PCFactorSetZeroPivot(PC pc,PetscReal zero)
 82: {

 88:   PetscTryMethod(pc,"PCFactorSetZeroPivot_C",(PC,PetscReal),(pc,zero));
 89:   return(0);
 90: }

 92: /*@
 93:    PCFactorSetShiftType - adds a particular type of quantity to the diagonal of the matrix during
 94:      numerical factorization, thus the matrix has nonzero pivots

 96:    Logically Collective on PC

 98:    Input Parameters:
 99: +  pc - the preconditioner context
100: -  shifttype - type of shift; one of MAT_SHIFT_NONE, MAT_SHIFT_NONZERO,  MAT_SHIFT_POSITIVE_DEFINITE, MAT_SHIFT_INBLOCKS

102:    Options Database Key:
103: .  -pc_factor_shift_type <shifttype> - Sets shift type or PETSC_DECIDE for the default; use '-help' for a list of available types

105:    Level: intermediate

107: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftAmount()
108: @*/
109: PetscErrorCode  PCFactorSetShiftType(PC pc,MatFactorShiftType shifttype)
110: {

116:   PetscTryMethod(pc,"PCFactorSetShiftType_C",(PC,MatFactorShiftType),(pc,shifttype));
117:   return(0);
118: }

120: /*@
121:    PCFactorSetShiftAmount - adds a quantity to the diagonal of the matrix during
122:      numerical factorization, thus the matrix has nonzero pivots

124:    Logically Collective on PC

126:    Input Parameters:
127: +  pc - the preconditioner context
128: -  shiftamount - amount of shift

130:    Options Database Key:
131: .  -pc_factor_shift_amount <shiftamount> - Sets shift amount or PETSC_DECIDE for the default

133:    Level: intermediate

135: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftType()
136: @*/
137: PetscErrorCode  PCFactorSetShiftAmount(PC pc,PetscReal shiftamount)
138: {

144:   PetscTryMethod(pc,"PCFactorSetShiftAmount_C",(PC,PetscReal),(pc,shiftamount));
145:   return(0);
146: }

148: /*
149:    PCFactorSetDropTolerance - The preconditioner will use an ILU
150:    based on a drop tolerance. (Under development)

152:    Logically Collective on PC

154:    Input Parameters:
155: +  pc - the preconditioner context
156: .  dt - the drop tolerance, try from 1.e-10 to .1
157: .  dtcol - tolerance for column pivot, good values [0.1 to 0.01]
158: -  maxrowcount - the max number of nonzeros allowed in a row, best value
159:                  depends on the number of nonzeros in row of original matrix

161:    Options Database Key:
162: .  -pc_factor_drop_tolerance <dt,dtcol,maxrowcount> - Sets drop tolerance

164:    Level: intermediate

166:       There are NO default values for the 3 parameters, you must set them with reasonable values for your
167:       matrix. We don't know how to compute reasonable values.

169: */
170: PetscErrorCode  PCFactorSetDropTolerance(PC pc,PetscReal dt,PetscReal dtcol,PetscInt maxrowcount)
171: {

178:   PetscTryMethod(pc,"PCFactorSetDropTolerance_C",(PC,PetscReal,PetscReal,PetscInt),(pc,dt,dtcol,maxrowcount));
179:   return(0);
180: }

182: /*@
183:    PCFactorGetZeroPivot - Gets the tolerance used to define a zero privot

185:    Not Collective

187:    Input Parameters:
188: .  pc - the preconditioner context

190:    Output Parameter:
191: .  pivot - the tolerance

193:    Level: intermediate


196: .seealso: PCFactorSetZeroPivot()
197: @*/
198: PetscErrorCode  PCFactorGetZeroPivot(PC pc,PetscReal *pivot)
199: {

204:   PetscUseMethod(pc,"PCFactorGetZeroPivot_C",(PC,PetscReal*),(pc,pivot));
205:   return(0);
206: }

208: /*@
209:    PCFactorGetShiftAmount - Gets the tolerance used to define a zero privot

211:    Not Collective

213:    Input Parameters:
214: .  pc - the preconditioner context

216:    Output Parameter:
217: .  shift - how much to shift the diagonal entry

219:    Level: intermediate


222: .seealso: PCFactorSetShiftAmount(), PCFactorSetShiftType(), PCFactorGetShiftType()
223: @*/
224: PetscErrorCode  PCFactorGetShiftAmount(PC pc,PetscReal *shift)
225: {

230:   PetscUseMethod(pc,"PCFactorGetShiftAmount_C",(PC,PetscReal*),(pc,shift));
231:   return(0);
232: }

234: /*@
235:    PCFactorGetShiftType - Gets the type of shift, if any, done when a zero pivot is detected

237:    Not Collective

239:    Input Parameters:
240: .  pc - the preconditioner context

242:    Output Parameter:
243: .  type - one of MAT_SHIFT_NONE, MAT_SHIFT_NONZERO,  MAT_SHIFT_POSITIVE_DEFINITE, or MAT_SHIFT_INBLOCKS

245:    Level: intermediate


248: .seealso: PCFactorSetShiftType(), MatFactorShiftType, PCFactorSetShiftAmount(), PCFactorGetShiftAmount()
249: @*/
250: PetscErrorCode  PCFactorGetShiftType(PC pc,MatFactorShiftType *type)
251: {

256:   PetscUseMethod(pc,"PCFactorGetShiftType_C",(PC,MatFactorShiftType*),(pc,type));
257:   return(0);
258: }

260: /*@
261:    PCFactorGetLevels - Gets the number of levels of fill to use.

263:    Logically Collective on PC

265:    Input Parameters:
266: .  pc - the preconditioner context

268:    Output Parameter:
269: .  levels - number of levels of fill

271:    Level: intermediate

273: @*/
274: PetscErrorCode  PCFactorGetLevels(PC pc,PetscInt *levels)
275: {

280:   PetscUseMethod(pc,"PCFactorGetLevels_C",(PC,PetscInt*),(pc,levels));
281:   return(0);
282: }

284: /*@
285:    PCFactorSetLevels - Sets the number of levels of fill to use.

287:    Logically Collective on PC

289:    Input Parameters:
290: +  pc - the preconditioner context
291: -  levels - number of levels of fill

293:    Options Database Key:
294: .  -pc_factor_levels <levels> - Sets fill level

296:    Level: intermediate

298: @*/
299: PetscErrorCode  PCFactorSetLevels(PC pc,PetscInt levels)
300: {

305:   if (levels < 0) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"negative levels");
307:   PetscTryMethod(pc,"PCFactorSetLevels_C",(PC,PetscInt),(pc,levels));
308:   return(0);
309: }

311: /*@
312:    PCFactorSetAllowDiagonalFill - Causes all diagonal matrix entries to be
313:    treated as level 0 fill even if there is no non-zero location.

315:    Logically Collective on PC

317:    Input Parameters:
318: +  pc - the preconditioner context
319: -  flg - PETSC_TRUE to turn on, PETSC_FALSE to turn off

321:    Options Database Key:
322: .  -pc_factor_diagonal_fill

324:    Notes:
325:    Does not apply with 0 fill.

327:    Level: intermediate

329: .seealso: PCFactorGetAllowDiagonalFill()

331: @*/
332: PetscErrorCode  PCFactorSetAllowDiagonalFill(PC pc,PetscBool flg)
333: {

338:   PetscTryMethod(pc,"PCFactorSetAllowDiagonalFill_C",(PC,PetscBool),(pc,flg));
339:   return(0);
340: }

342: /*@
343:    PCFactorGetAllowDiagonalFill - Determines if all diagonal matrix entries are
344:        treated as level 0 fill even if there is no non-zero location.

346:    Logically Collective on PC

348:    Input Parameter:
349: .  pc - the preconditioner context

351:    Output Parameter:
352: .   flg - PETSC_TRUE to turn on, PETSC_FALSE to turn off

354:    Options Database Key:
355: .  -pc_factor_diagonal_fill

357:    Notes:
358:    Does not apply with 0 fill.

360:    Level: intermediate

362: .seealso: PCFactorSetAllowDiagonalFill()

364: @*/
365: PetscErrorCode  PCFactorGetAllowDiagonalFill(PC pc,PetscBool *flg)
366: {

371:   PetscUseMethod(pc,"PCFactorGetAllowDiagonalFill_C",(PC,PetscBool*),(pc,flg));
372:   return(0);
373: }

375: /*@
376:    PCFactorReorderForNonzeroDiagonal - reorders rows/columns of matrix to remove zeros from diagonal

378:    Logically Collective on PC

380:    Input Parameters:
381: +  pc - the preconditioner context
382: -  tol - diagonal entries smaller than this in absolute value are considered zero

384:    Options Database Key:
385: .  -pc_factor_nonzeros_along_diagonal <tol>

387:    Level: intermediate

389: .seealso: PCFactorSetFill(), PCFactorSetShiftNonzero(), PCFactorSetZeroPivot(), MatReorderForNonzeroDiagonal()
390: @*/
391: PetscErrorCode  PCFactorReorderForNonzeroDiagonal(PC pc,PetscReal rtol)
392: {

398:   PetscTryMethod(pc,"PCFactorReorderForNonzeroDiagonal_C",(PC,PetscReal),(pc,rtol));
399:   return(0);
400: }

402: /*@C
403:    PCFactorSetMatSolverType - sets the software that is used to perform the factorization

405:    Logically Collective on PC

407:    Input Parameters:
408: +  pc - the preconditioner context
409: -  stype - for example, superlu, superlu_dist

411:    Options Database Key:
412: .  -pc_factor_mat_solver_type <stype> - petsc, superlu, superlu_dist, mumps, cusparse

414:    Level: intermediate

416:    Note:
417:      By default this will use the PETSc factorization if it exists


420: .seealso: MatGetFactor(), MatSolverType, PCFactorGetMatSolverType()

422: @*/
423: PetscErrorCode  PCFactorSetMatSolverType(PC pc,MatSolverType stype)
424: {

429:   PetscTryMethod(pc,"PCFactorSetMatSolverType_C",(PC,MatSolverType),(pc,stype));
430:   return(0);
431: }

433: /*@C
434:    PCFactorGetMatSolverType - gets the software that is used to perform the factorization

436:    Not Collective

438:    Input Parameter:
439: .  pc - the preconditioner context

441:    Output Parameter:
442: .   stype - for example, superlu, superlu_dist (NULL if the PC does not have a solver package)

444:    Level: intermediate


447: .seealso: MatGetFactor(), MatSolverType, PCFactorGetMatSolverType()

449: @*/
450: PetscErrorCode  PCFactorGetMatSolverType(PC pc,MatSolverType *stype)
451: {
452:   PetscErrorCode ierr,(*f)(PC,MatSolverType*);

456:   PetscObjectQueryFunction((PetscObject)pc,"PCFactorGetMatSolverType_C",&f);
457:   if (f) {
458:     (*f)(pc,stype);
459:   } else {
460:     *stype = NULL;
461:   }
462:   return(0);
463: }

465: /*@
466:    PCFactorSetFill - Indicate the amount of fill you expect in the factored matrix,
467:    fill = number nonzeros in factor/number nonzeros in original matrix.

469:    Not Collective, each process can expect a different amount of fill

471:    Input Parameters:
472: +  pc - the preconditioner context
473: -  fill - amount of expected fill

475:    Options Database Key:
476: .  -pc_factor_fill <fill> - Sets fill amount

478:    Level: intermediate

480:    Note:
481:    For sparse matrix factorizations it is difficult to predict how much
482:    fill to expect. By running with the option -info PETSc will print the
483:    actual amount of fill used; allowing you to set the value accurately for
484:    future runs. Default PETSc uses a value of 5.0

486:    This parameter has NOTHING to do with the levels-of-fill of ILU(). That is set with PCFactorSetLevels() or -pc_factor_levels.


489: @*/
490: PetscErrorCode  PCFactorSetFill(PC pc,PetscReal fill)
491: {

496:   if (fill < 1.0) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Fill factor cannot be less then 1.0");
497:   PetscTryMethod(pc,"PCFactorSetFill_C",(PC,PetscReal),(pc,fill));
498:   return(0);
499: }

501: /*@
502:    PCFactorSetUseInPlace - Tells the system to do an in-place factorization.
503:    For dense matrices, this enables the solution of much larger problems.
504:    For sparse matrices the factorization cannot be done truly in-place
505:    so this does not save memory during the factorization, but after the matrix
506:    is factored, the original unfactored matrix is freed, thus recovering that
507:    space.

509:    Logically Collective on PC

511:    Input Parameters:
512: +  pc - the preconditioner context
513: -  flg - PETSC_TRUE to enable, PETSC_FALSE to disable

515:    Options Database Key:
516: .  -pc_factor_in_place <true,false>- Activate/deactivate in-place factorization

518:    Notes:
519:    PCFactorSetUseInplace() can only be used with the KSP method KSPPREONLY or when
520:    a different matrix is provided for the multiply and the preconditioner in
521:    a call to KSPSetOperators().
522:    This is because the Krylov space methods require an application of the
523:    matrix multiplication, which is not possible here because the matrix has
524:    been factored in-place, replacing the original matrix.

526:    Level: intermediate

528: .seealso: PCFactorGetUseInPlace()
529: @*/
530: PetscErrorCode  PCFactorSetUseInPlace(PC pc,PetscBool flg)
531: {

536:   PetscTryMethod(pc,"PCFactorSetUseInPlace_C",(PC,PetscBool),(pc,flg));
537:   return(0);
538: }

540: /*@
541:    PCFactorGetUseInPlace - Determines if an in-place factorization is being used.

543:    Logically Collective on PC

545:    Input Parameter:
546: .  pc - the preconditioner context

548:    Output Parameter:
549: .  flg - PETSC_TRUE to enable, PETSC_FALSE to disable

551:    Level: intermediate

553: .seealso: PCFactorSetUseInPlace()
554: @*/
555: PetscErrorCode  PCFactorGetUseInPlace(PC pc,PetscBool *flg)
556: {

561:   PetscUseMethod(pc,"PCFactorGetUseInPlace_C",(PC,PetscBool*),(pc,flg));
562:   return(0);
563: }

565: /*@C
566:     PCFactorSetMatOrderingType - Sets the ordering routine (to reduce fill) to
567:     be used in the LU factorization.

569:     Logically Collective on PC

571:     Input Parameters:
572: +   pc - the preconditioner context
573: -   ordering - the matrix ordering name, for example, MATORDERINGND or MATORDERINGRCM

575:     Options Database Key:
576: .   -pc_factor_mat_ordering_type <nd,rcm,...> - Sets ordering routine

578:     Level: intermediate

580:     Notes:
581:     nested dissection is used by default

583:     For Cholesky and ICC and the SBAIJ format reorderings are not available,
584:     since only the upper triangular part of the matrix is stored. You can use the
585:     SeqAIJ format in this case to get reorderings.

587: @*/
588: PetscErrorCode  PCFactorSetMatOrderingType(PC pc,MatOrderingType ordering)
589: {

594:   PetscTryMethod(pc,"PCFactorSetMatOrderingType_C",(PC,MatOrderingType),(pc,ordering));
595:   return(0);
596: }

598: /*@
599:     PCFactorSetColumnPivot - Determines when column pivoting is done during matrix factorization.
600:       For PETSc dense matrices column pivoting is always done, for PETSc sparse matrices
601:       it is never done. For the MATLAB and SuperLU factorization this is used.

603:     Logically Collective on PC

605:     Input Parameters:
606: +   pc - the preconditioner context
607: -   dtcol - 0.0 implies no pivoting, 1.0 complete pivoting (slower, requires more memory but more stable)

609:     Options Database Key:
610: .   -pc_factor_pivoting <dtcol>

612:     Level: intermediate

614: .seealso: PCILUSetMatOrdering(), PCFactorSetPivotInBlocks()
615: @*/
616: PetscErrorCode  PCFactorSetColumnPivot(PC pc,PetscReal dtcol)
617: {

623:   PetscTryMethod(pc,"PCFactorSetColumnPivot_C",(PC,PetscReal),(pc,dtcol));
624:   return(0);
625: }

627: /*@
628:     PCFactorSetPivotInBlocks - Determines if pivoting is done while factoring each block
629:       with BAIJ or SBAIJ matrices

631:     Logically Collective on PC

633:     Input Parameters:
634: +   pc - the preconditioner context
635: -   pivot - PETSC_TRUE or PETSC_FALSE

637:     Options Database Key:
638: .   -pc_factor_pivot_in_blocks <true,false>

640:     Level: intermediate

642: .seealso: PCILUSetMatOrdering(), PCFactorSetColumnPivot()
643: @*/
644: PetscErrorCode  PCFactorSetPivotInBlocks(PC pc,PetscBool pivot)
645: {

651:   PetscTryMethod(pc,"PCFactorSetPivotInBlocks_C",(PC,PetscBool),(pc,pivot));
652:   return(0);
653: }

655: /*@
656:    PCFactorSetReuseFill - When matrices with different nonzero structure are factored,
657:    this causes later ones to use the fill ratio computed in the initial factorization.

659:    Logically Collective on PC

661:    Input Parameters:
662: +  pc - the preconditioner context
663: -  flag - PETSC_TRUE to reuse else PETSC_FALSE

665:    Options Database Key:
666: .  -pc_factor_reuse_fill - Activates PCFactorSetReuseFill()

668:    Level: intermediate

670: .seealso: PCFactorSetReuseOrdering()
671: @*/
672: PetscErrorCode  PCFactorSetReuseFill(PC pc,PetscBool flag)
673: {

679:   PetscTryMethod(pc,"PCFactorSetReuseFill_C",(PC,PetscBool),(pc,flag));
680:   return(0);
681: }

683: PetscErrorCode PCFactorInitialize(PC pc)
684: {
686:   PC_Factor       *fact = (PC_Factor*)pc->data;

689:   MatFactorInfoInitialize(&fact->info);
690:   fact->info.shifttype       = (PetscReal)MAT_SHIFT_NONE;
691:   fact->info.shiftamount     = 100.0*PETSC_MACHINE_EPSILON;
692:   fact->info.zeropivot       = 100.0*PETSC_MACHINE_EPSILON;
693:   fact->info.pivotinblocks   = 1.0;
694:   pc->ops->getfactoredmatrix = PCFactorGetMatrix_Factor;

696:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetZeroPivot_C",PCFactorSetZeroPivot_Factor);
697:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetZeroPivot_C",PCFactorGetZeroPivot_Factor);
698:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftType_C",PCFactorSetShiftType_Factor);
699:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetShiftType_C",PCFactorGetShiftType_Factor);
700:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftAmount_C",PCFactorSetShiftAmount_Factor);
701:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetShiftAmount_C",PCFactorGetShiftAmount_Factor);
702:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetMatSolverType_C",PCFactorGetMatSolverType_Factor);
703:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatSolverType_C",PCFactorSetMatSolverType_Factor);
704:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUpMatSolverType_C",PCFactorSetUpMatSolverType_Factor);
705:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetFill_C",PCFactorSetFill_Factor);
706:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatOrderingType_C",PCFactorSetMatOrderingType_Factor);
707:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetLevels_C",PCFactorSetLevels_Factor);
708:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetLevels_C",PCFactorGetLevels_Factor);
709:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetAllowDiagonalFill_C",PCFactorSetAllowDiagonalFill_Factor);
710:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetAllowDiagonalFill_C",PCFactorGetAllowDiagonalFill_Factor);
711:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetPivotInBlocks_C",PCFactorSetPivotInBlocks_Factor);
712:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUseInPlace_C",PCFactorSetUseInPlace_Factor);
713:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetUseInPlace_C",PCFactorGetUseInPlace_Factor);
714:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseOrdering_C",PCFactorSetReuseOrdering_Factor);
715:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseFill_C",PCFactorSetReuseFill_Factor);
716:   return(0);
717: }