Actual source code: factor.c

  1: #include <../src/ksp/pc/impls/factor/factor.h>
  2: #include <petsc/private/matimpl.h>

  4: /*
  5:     If an ordering is not yet set and the matrix is available determine a default ordering
  6: */
  7: PetscErrorCode PCFactorSetDefaultOrdering_Factor(PC pc)
  8: {
  9:   PetscBool   foundmtype, flg, destroy = PETSC_FALSE;
 10:   const char *prefix;

 12:   PetscFunctionBegin;
 13:   if (pc->pmat) {
 14:     PetscCall(PCGetOptionsPrefix(pc, &prefix));
 15:     PetscCall(MatSetOptionsPrefixFactor(pc->pmat, prefix));
 16:     PC_Factor *fact = (PC_Factor *)pc->data;
 17:     PetscCall(MatSolverTypeGet(fact->solvertype, ((PetscObject)pc->pmat)->type_name, fact->factortype, NULL, &foundmtype, NULL));
 18:     if (foundmtype) {
 19:       if (!fact->fact) {
 20:         /* factored matrix is not present at this point, we want to create it during PCSetUp.
 21:            Since this can be called from setfromoptions, we destroy it when we are done with it */
 22:         PetscCall(MatGetFactor(pc->pmat, fact->solvertype, fact->factortype, &fact->fact));
 23:         destroy = PETSC_TRUE;
 24:       }
 25:       if (!fact->fact) PetscFunctionReturn(PETSC_SUCCESS);
 26:       if (!fact->fact->assembled) {
 27:         PetscCall(PetscStrcmp(fact->solvertype, fact->fact->solvertype, &flg));
 28:         if (!flg) {
 29:           Mat B;
 30:           PetscCall(MatGetFactor(pc->pmat, fact->solvertype, fact->factortype, &B));
 31:           PetscCall(MatHeaderReplace(fact->fact, &B));
 32:         }
 33:       }
 34:       if (!fact->ordering) {
 35:         PetscBool       canuseordering;
 36:         MatOrderingType otype;

 38:         PetscCall(MatFactorGetCanUseOrdering(fact->fact, &canuseordering));
 39:         if (canuseordering) {
 40:           PetscCall(MatFactorGetPreferredOrdering(fact->fact, fact->factortype, &otype));
 41:         } else otype = MATORDERINGEXTERNAL;
 42:         PetscCall(PetscStrallocpy(otype, (char **)&fact->ordering));
 43:       }
 44:       if (destroy) PetscCall(MatDestroy(&fact->fact));
 45:     }
 46:   }
 47:   PetscFunctionReturn(PETSC_SUCCESS);
 48: }

 50: static PetscErrorCode PCFactorSetReuseOrdering_Factor(PC pc, PetscBool flag)
 51: {
 52:   PC_Factor *lu = (PC_Factor *)pc->data;

 54:   PetscFunctionBegin;
 55:   lu->reuseordering = flag;
 56:   PetscFunctionReturn(PETSC_SUCCESS);
 57: }

 59: static PetscErrorCode PCFactorSetReuseFill_Factor(PC pc, PetscBool flag)
 60: {
 61:   PC_Factor *lu = (PC_Factor *)pc->data;

 63:   PetscFunctionBegin;
 64:   lu->reusefill = flag;
 65:   PetscFunctionReturn(PETSC_SUCCESS);
 66: }

 68: static PetscErrorCode PCFactorSetUseInPlace_Factor(PC pc, PetscBool flg)
 69: {
 70:   PC_Factor *dir = (PC_Factor *)pc->data;

 72:   PetscFunctionBegin;
 73:   dir->inplace = flg;
 74:   PetscFunctionReturn(PETSC_SUCCESS);
 75: }

 77: static PetscErrorCode PCFactorGetUseInPlace_Factor(PC pc, PetscBool *flg)
 78: {
 79:   PC_Factor *dir = (PC_Factor *)pc->data;

 81:   PetscFunctionBegin;
 82:   *flg = dir->inplace;
 83:   PetscFunctionReturn(PETSC_SUCCESS);
 84: }

 86: /*@
 87:   PCFactorSetUpMatSolverType - Can be called after `KSPSetOperators()` or `PCSetOperators()`, causes `MatGetFactor()` to be called so then one may
 88:   set the options for that particular factorization object.

 90:   Input Parameter:
 91: . pc - the preconditioner context

 93:   Note:
 94:   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.
 95:   This function raises an error if the requested combination of solver package and matrix type is not supported.

 97:   Level: intermediate

 99: .seealso: [](ch_ksp), `PCCHOLESKY`, `PCLU`, `PCFactorSetMatSolverType()`, `PCFactorGetMatrix()`
100: @*/
101: PetscErrorCode PCFactorSetUpMatSolverType(PC pc)
102: {
103:   PetscFunctionBegin;
105:   PetscTryMethod(pc, "PCFactorSetUpMatSolverType_C", (PC), (pc));
106:   PetscFunctionReturn(PETSC_SUCCESS);
107: }

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

112:   Logically Collective

114:   Input Parameters:
115: + pc   - the preconditioner context
116: - zero - all pivots smaller than this will be considered zero

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

121:   Level: intermediate

123: .seealso: [](ch_ksp), `PCCHOLESKY`, `PCLU`, `PCFactorSetShiftType()`, `PCFactorSetShiftAmount()`
124: @*/
125: PetscErrorCode PCFactorSetZeroPivot(PC pc, PetscReal zero)
126: {
127:   PetscFunctionBegin;
130:   PetscTryMethod(pc, "PCFactorSetZeroPivot_C", (PC, PetscReal), (pc, zero));
131:   PetscFunctionReturn(PETSC_SUCCESS);
132: }

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

138:   Logically Collective

140:   Input Parameters:
141: + pc        - the preconditioner context
142: - shifttype - type of shift; one of `MAT_SHIFT_NONE`, `MAT_SHIFT_NONZERO`, `MAT_SHIFT_POSITIVE_DEFINITE`, `MAT_SHIFT_INBLOCKS`

144:   Options Database Key:
145: . -pc_factor_shift_type <shifttype> - Sets shift type; use '-help' for a list of available types

147:   Level: intermediate

149: .seealso: [](ch_ksp), `PCCHOLESKY`, `PCLU`, `PCFactorSetZeroPivot()`, `PCFactorSetShiftAmount()`
150: @*/
151: PetscErrorCode PCFactorSetShiftType(PC pc, MatFactorShiftType shifttype)
152: {
153:   PetscFunctionBegin;
156:   PetscTryMethod(pc, "PCFactorSetShiftType_C", (PC, MatFactorShiftType), (pc, shifttype));
157:   PetscFunctionReturn(PETSC_SUCCESS);
158: }

160: /*@
161:   PCFactorSetShiftAmount - adds a quantity to the diagonal of the matrix during
162:   numerical factorization, thus the matrix has nonzero pivots

164:   Logically Collective

166:   Input Parameters:
167: + pc          - the preconditioner context
168: - shiftamount - amount of shift or `PETSC_DECIDE` for the default

170:   Options Database Key:
171: . -pc_factor_shift_amount <shiftamount> - Sets shift amount or -1 for the default

173:   Level: intermediate

175: .seealso: [](ch_ksp), `PCCHOLESKY`, `PCLU`, ``PCFactorSetZeroPivot()`, `PCFactorSetShiftType()`
176: @*/
177: PetscErrorCode PCFactorSetShiftAmount(PC pc, PetscReal shiftamount)
178: {
179:   PetscFunctionBegin;
182:   PetscTryMethod(pc, "PCFactorSetShiftAmount_C", (PC, PetscReal), (pc, shiftamount));
183:   PetscFunctionReturn(PETSC_SUCCESS);
184: }

186: /*@
187:   PCFactorSetDropTolerance - The preconditioner will use an `PCILU`
188:   based on a drop tolerance.

190:   Logically Collective

192:   Input Parameters:
193: + pc          - the preconditioner context
194: . dt          - the drop tolerance, try from 1.e-10 to .1
195: . dtcol       - tolerance for column pivot, good values [0.1 to 0.01]
196: - maxrowcount - the max number of nonzeros allowed in a row, best value
197:                  depends on the number of nonzeros in row of original matrix

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

202:   Level: intermediate

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

208: .seealso: [](ch_ksp), `PCILU`
209: @*/
210: PetscErrorCode PCFactorSetDropTolerance(PC pc, PetscReal dt, PetscReal dtcol, PetscInt maxrowcount)
211: {
212:   PetscFunctionBegin;
216:   PetscTryMethod(pc, "PCFactorSetDropTolerance_C", (PC, PetscReal, PetscReal, PetscInt), (pc, dt, dtcol, maxrowcount));
217:   PetscFunctionReturn(PETSC_SUCCESS);
218: }

220: /*@
221:   PCFactorGetZeroPivot - Gets the tolerance used to define a zero privot

223:   Not Collective

225:   Input Parameter:
226: . pc - the preconditioner context

228:   Output Parameter:
229: . pivot - the tolerance

231:   Level: intermediate

233: .seealso: [](ch_ksp), `PCLU`, `PCCHOLESKY`, `PCFactorSetZeroPivot()`
234: @*/
235: PetscErrorCode PCFactorGetZeroPivot(PC pc, PetscReal *pivot)
236: {
237:   PetscFunctionBegin;
239:   PetscUseMethod(pc, "PCFactorGetZeroPivot_C", (PC, PetscReal *), (pc, pivot));
240:   PetscFunctionReturn(PETSC_SUCCESS);
241: }

243: /*@
244:   PCFactorGetShiftAmount - Gets the tolerance used to define a zero privot

246:   Not Collective

248:   Input Parameter:
249: . pc - the preconditioner context

251:   Output Parameter:
252: . shift - how much to shift the diagonal entry

254:   Level: intermediate

256: .seealso: [](ch_ksp), `PCLU`, `PCCHOLESKY`, `PCFactorSetShiftAmount()`, `PCFactorSetShiftType()`, `PCFactorGetShiftType()`
257: @*/
258: PetscErrorCode PCFactorGetShiftAmount(PC pc, PetscReal *shift)
259: {
260:   PetscFunctionBegin;
262:   PetscUseMethod(pc, "PCFactorGetShiftAmount_C", (PC, PetscReal *), (pc, shift));
263:   PetscFunctionReturn(PETSC_SUCCESS);
264: }

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

269:   Not Collective

271:   Input Parameter:
272: . pc - the preconditioner context

274:   Output Parameter:
275: . type - one of `MAT_SHIFT_NONE`, `MAT_SHIFT_NONZERO`, `MAT_SHIFT_POSITIVE_DEFINITE`, or `MAT_SHIFT_INBLOCKS`

277:   Level: intermediate

279: .seealso: [](ch_ksp), `PCLU`, `PCCHOLESKY`, `PCFactorSetShiftType()`, `MatFactorShiftType`, `PCFactorSetShiftAmount()`, `PCFactorGetShiftAmount()`
280: @*/
281: PetscErrorCode PCFactorGetShiftType(PC pc, MatFactorShiftType *type)
282: {
283:   PetscFunctionBegin;
285:   PetscUseMethod(pc, "PCFactorGetShiftType_C", (PC, MatFactorShiftType *), (pc, type));
286:   PetscFunctionReturn(PETSC_SUCCESS);
287: }

289: /*@
290:   PCFactorGetLevels - Gets the number of levels of fill to use.

292:   Logically Collective

294:   Input Parameter:
295: . pc - the preconditioner context

297:   Output Parameter:
298: . levels - number of levels of fill

300:   Level: intermediate

302: .seealso: [](ch_ksp), `PCILU`, `PCICC`, `PCFactorSetLevels()`
303: @*/
304: PetscErrorCode PCFactorGetLevels(PC pc, PetscInt *levels)
305: {
306:   PetscFunctionBegin;
308:   PetscUseMethod(pc, "PCFactorGetLevels_C", (PC, PetscInt *), (pc, levels));
309:   PetscFunctionReturn(PETSC_SUCCESS);
310: }

312: /*@
313:   PCFactorSetLevels - Sets the number of levels of fill to use.

315:   Logically Collective

317:   Input Parameters:
318: + pc     - the preconditioner context
319: - levels - number of levels of fill

321:   Options Database Key:
322: . -pc_factor_levels <levels> - Sets fill level

324:   Level: intermediate

326: .seealso: [](ch_ksp), `PCILU`, `PCICC`, `PCFactorGetLevels()`
327: @*/
328: PetscErrorCode PCFactorSetLevels(PC pc, PetscInt levels)
329: {
330:   PetscFunctionBegin;
332:   PetscCheck(levels >= 0, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "negative levels");
334:   PetscTryMethod(pc, "PCFactorSetLevels_C", (PC, PetscInt), (pc, levels));
335:   PetscFunctionReturn(PETSC_SUCCESS);
336: }

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

342:   Logically Collective

344:   Input Parameters:
345: + pc  - the preconditioner context
346: - flg - `PETSC_TRUE` to turn on, `PETSC_FALSE` to turn off

348:   Options Database Key:
349: . -pc_factor_diagonal_fill <bool> - allow the diagonal fill

351:   Note:
352:   Does not apply with 0 fill.

354:   Level: intermediate

356: .seealso: [](ch_ksp), `PCILU`, `PCICC`, `PCFactorGetAllowDiagonalFill()`
357: @*/
358: PetscErrorCode PCFactorSetAllowDiagonalFill(PC pc, PetscBool flg)
359: {
360:   PetscFunctionBegin;
362:   PetscTryMethod(pc, "PCFactorSetAllowDiagonalFill_C", (PC, PetscBool), (pc, flg));
363:   PetscFunctionReturn(PETSC_SUCCESS);
364: }

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

370:   Logically Collective

372:   Input Parameter:
373: . pc - the preconditioner context

375:   Output Parameter:
376: . flg - `PETSC_TRUE` to turn on, `PETSC_FALSE` to turn off

378:   Note:
379:   Does not apply with 0 fill.

381:   Level: intermediate

383: .seealso: [](ch_ksp), `PCILU`, `PCICC`, `PCFactorSetAllowDiagonalFill()`
384: @*/
385: PetscErrorCode PCFactorGetAllowDiagonalFill(PC pc, PetscBool *flg)
386: {
387:   PetscFunctionBegin;
389:   PetscUseMethod(pc, "PCFactorGetAllowDiagonalFill_C", (PC, PetscBool *), (pc, flg));
390:   PetscFunctionReturn(PETSC_SUCCESS);
391: }

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

396:   Logically Collective

398:   Input Parameters:
399: + pc   - the preconditioner context
400: - rtol - diagonal entries smaller than this in absolute value are considered zero

402:   Options Database Key:
403: . -pc_factor_nonzeros_along_diagonal <tol> - perform the reordering with the given tolerance

405:   Level: intermediate

407: .seealso: [](ch_ksp), `PCILU`, `PCICC`, `PCFactorSetFill()`, `PCFactorSetShiftNonzero()`, `PCFactorSetZeroPivot()`, `MatReorderForNonzeroDiagonal()`
408: @*/
409: PetscErrorCode PCFactorReorderForNonzeroDiagonal(PC pc, PetscReal rtol)
410: {
411:   PetscFunctionBegin;
414:   PetscTryMethod(pc, "PCFactorReorderForNonzeroDiagonal_C", (PC, PetscReal), (pc, rtol));
415:   PetscFunctionReturn(PETSC_SUCCESS);
416: }

418: /*@C
419:   PCFactorSetMatSolverType - sets the solver package that is used to perform the factorization

421:   Logically Collective

423:   Input Parameters:
424: + pc    - the preconditioner context
425: - stype - for example, `MATSOLVERSUPERLU`, `MATSOLVERSUPERLU_DIST`, `MATSOLVERMUMPS`

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

430:   Level: intermediate

432:   Note:
433:   The default type is set by searching for available types based on the order of the calls to `MatSolverTypeRegister()` in `MatInitializePackage()`.
434:   Since different PETSc configurations may have different external solvers, seemingly identical runs with different PETSc configurations may use a different solver.
435:   For example if one configuration had --download-mumps while a different one had --download-superlu_dist.

437: .seealso: [](ch_ksp), `PCLU`, `PCCHOLESKY`, `MatGetFactor()`, `MatSolverType`, `PCFactorGetMatSolverType()`, `MatSolverTypeRegister()`,
438:           `MatInitializePackage()`, `MATSOLVERSUPERLU`, `MATSOLVERSUPERLU_DIST`, `MATSOLVERMUMPS`, `MatSolverTypeGet()`
439: @*/
440: PetscErrorCode PCFactorSetMatSolverType(PC pc, MatSolverType stype)
441: {
442:   PetscFunctionBegin;
444:   PetscTryMethod(pc, "PCFactorSetMatSolverType_C", (PC, MatSolverType), (pc, stype));
445:   PetscFunctionReturn(PETSC_SUCCESS);
446: }

448: /*@C
449:   PCFactorGetMatSolverType - gets the solver package that is used to perform the factorization

451:   Not Collective

453:   Input Parameter:
454: . pc - the preconditioner context

456:   Output Parameter:
457: . stype - for example, `MATSOLVERSUPERLU`, `MATSOLVERSUPERLU_DIST`, `MATSOLVERMUMPS`

459:   Level: intermediate

461: .seealso: [](ch_ksp), `PCLU`, `PCCHOLESKY`, `MatGetFactor()`, `MatSolverType`, `MATSOLVERSUPERLU`,
462: `MATSOLVERSUPERLU_DIST`, `MATSOLVERMUMPS`
463: @*/
464: PetscErrorCode PCFactorGetMatSolverType(PC pc, MatSolverType *stype)
465: {
466:   PetscErrorCode (*f)(PC, MatSolverType *);

468:   PetscFunctionBegin;
470:   PetscAssertPointer(stype, 2);
471:   PetscCall(PetscObjectQueryFunction((PetscObject)pc, "PCFactorGetMatSolverType_C", &f));
472:   if (f) PetscCall((*f)(pc, stype));
473:   else *stype = NULL;
474:   PetscFunctionReturn(PETSC_SUCCESS);
475: }

477: /*@
478:   PCFactorSetFill - Indicate the amount of fill you expect in the factored matrix,
479:   fill = number nonzeros in factor/number nonzeros in original matrix.

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

483:   Input Parameters:
484: + pc   - the preconditioner context
485: - fill - amount of expected fill

487:   Options Database Key:
488: . -pc_factor_fill <fill> - Sets fill amount

490:   Level: intermediate

492:   Notes:
493:   For sparse matrix factorizations it is difficult to predict how much
494:   fill to expect. By running with the option -info PETSc will print the
495:   actual amount of fill used; allowing you to set the value accurately for
496:   future runs. Default PETSc uses a value of 5.0

498:   This is ignored for most solver packages

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

502: .seealso: [](ch_ksp), `PCLU`, `PCCHOLESKY`, `PCILU`, `PCICC`, `PCFactorSetReuseFill()`
503: @*/
504: PetscErrorCode PCFactorSetFill(PC pc, PetscReal fill)
505: {
506:   PetscFunctionBegin;
508:   PetscCheck(fill >= 1.0, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Fill factor cannot be less than 1.0");
509:   PetscTryMethod(pc, "PCFactorSetFill_C", (PC, PetscReal), (pc, fill));
510:   PetscFunctionReturn(PETSC_SUCCESS);
511: }

513: /*@
514:   PCFactorSetUseInPlace - Tells the preconditioner to do an in-place factorization.

516:   Logically Collective

518:   Input Parameters:
519: + pc  - the preconditioner context
520: - flg - `PETSC_TRUE` to enable, `PETSC_FALSE` to disable

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

525:   Note:
526:   For dense matrices, this enables the solution of much larger problems.
527:   For sparse matrices the factorization cannot be done truly in-place
528:   so this does not save memory during the factorization, but after the matrix
529:   is factored, the original unfactored matrix is freed, thus recovering that
530:   space. For ICC(0) and ILU(0) with the default natural ordering the factorization is done efficiently in-place.

532:   `PCFactorSetUseInplace()` can only be used with the `KSP` method `KSPPREONLY` or when
533:   a different matrix is provided for the multiply and the preconditioner in
534:   a call to `KSPSetOperators()`.
535:   This is because the Krylov space methods require an application of the
536:   matrix multiplication, which is not possible here because the matrix has
537:   been factored in-place, replacing the original matrix.

539:   Level: intermediate

541: .seealso: [](ch_ksp), `PC`, `Mat`, `PCLU`, `PCCHOLESKY`, `PCILU`, `PCICC`, `PCFactorGetUseInPlace()`
542: @*/
543: PetscErrorCode PCFactorSetUseInPlace(PC pc, PetscBool flg)
544: {
545:   PetscFunctionBegin;
547:   PetscTryMethod(pc, "PCFactorSetUseInPlace_C", (PC, PetscBool), (pc, flg));
548:   PetscFunctionReturn(PETSC_SUCCESS);
549: }

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

554:   Logically Collective

556:   Input Parameter:
557: . pc - the preconditioner context

559:   Output Parameter:
560: . flg - `PETSC_TRUE` to enable, `PETSC_FALSE` to disable

562:   Level: intermediate

564: .seealso: [](ch_ksp), `PCLU`, `PCCHOLESKY`, `PCILU`, `PCICC`, `PCFactorSetUseInPlace()`
565: @*/
566: PetscErrorCode PCFactorGetUseInPlace(PC pc, PetscBool *flg)
567: {
568:   PetscFunctionBegin;
570:   PetscUseMethod(pc, "PCFactorGetUseInPlace_C", (PC, PetscBool *), (pc, flg));
571:   PetscFunctionReturn(PETSC_SUCCESS);
572: }

574: /*@C
575:   PCFactorSetMatOrderingType - Sets the ordering routine (to reduce fill) to
576:   be used in the `PCLU`, `PCCHOLESKY`, `PCILU`,  or `PCICC` preconditioners

578:   Logically Collective

580:   Input Parameters:
581: + pc       - the preconditioner context
582: - ordering - the matrix ordering name, for example, `MATORDERINGND` or `MATORDERINGRCM`

584:   Options Database Key:
585: . -pc_factor_mat_ordering_type <nd,rcm,...,external> - Sets ordering routine

587:   Level: intermediate

589:   Notes:
590:   Nested dissection is used by default for some of PETSc's sparse matrix formats

592:   For `PCCHOLESKY` and `PCICC` and the `MATSBAIJ` format the only reordering available is natural since only the upper half of the matrix is stored
593:   and reordering this matrix is very expensive.

595:   You can use a `MATSEQAIJ` matrix with Cholesky and ICC and use any ordering.

597:   `MATORDERINGEXTERNAL` means PETSc will not compute an ordering and the package will use its own ordering, usable with `MATSOLVERCHOLMOD`, `MATSOLVERUMFPACK`, and others.

599: .seealso: [](ch_ksp), `PCLU`, `PCCHOLESKY`, `PCILU`, `PCICC`, `MatOrderingType`, `MATORDERINGEXTERNAL`, `MATORDERINGND`, `MATORDERINGRCM`
600: @*/
601: PetscErrorCode PCFactorSetMatOrderingType(PC pc, MatOrderingType ordering)
602: {
603:   PetscFunctionBegin;
605:   PetscTryMethod(pc, "PCFactorSetMatOrderingType_C", (PC, MatOrderingType), (pc, ordering));
606:   PetscFunctionReturn(PETSC_SUCCESS);
607: }

609: /*@
610:   PCFactorSetColumnPivot - Determines when column pivoting is done during matrix factorization.
611:   For PETSc dense matrices column pivoting is always done, for PETSc sparse matrices
612:   it is never done. For the MATLAB and `MATSOLVERSUPERLU` factorization this is used.

614:   Logically Collective

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

620:   Options Database Key:
621: . -pc_factor_pivoting <dtcol> - perform the pivoting with the given tolerance

623:   Level: intermediate

625: .seealso: [](ch_ksp), `PCLU`, `PCCHOLESKY`, `PCILU`, `PCICC`, `PCILUSetMatOrdering()`, `PCFactorSetPivotInBlocks()`
626: @*/
627: PetscErrorCode PCFactorSetColumnPivot(PC pc, PetscReal dtcol)
628: {
629:   PetscFunctionBegin;
632:   PetscTryMethod(pc, "PCFactorSetColumnPivot_C", (PC, PetscReal), (pc, dtcol));
633:   PetscFunctionReturn(PETSC_SUCCESS);
634: }

636: /*@
637:   PCFactorSetPivotInBlocks - Determines if pivoting is done while factoring each block
638:   with `MATBAIJ` or `MATSBAIJ` matrices

640:   Logically Collective

642:   Input Parameters:
643: + pc    - the preconditioner context
644: - pivot - `PETSC_TRUE` or `PETSC_FALSE`

646:   Options Database Key:
647: . -pc_factor_pivot_in_blocks <true,false> - Pivot inside matrix dense blocks for `MATBAIJ` and `MATSBAIJ`

649:   Level: intermediate

651: .seealso: [](ch_ksp), `PCLU`, `PCCHOLESKY`, `PCILU`, `PCICC`, `PCILUSetMatOrdering()`, `PCFactorSetColumnPivot()`
652: @*/
653: PetscErrorCode PCFactorSetPivotInBlocks(PC pc, PetscBool pivot)
654: {
655:   PetscFunctionBegin;
658:   PetscTryMethod(pc, "PCFactorSetPivotInBlocks_C", (PC, PetscBool), (pc, pivot));
659:   PetscFunctionReturn(PETSC_SUCCESS);
660: }

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

666:   Logically Collective

668:   Input Parameters:
669: + pc   - the preconditioner context
670: - flag - `PETSC_TRUE` to reuse else `PETSC_FALSE`

672:   Options Database Key:
673: . -pc_factor_reuse_fill - Activates `PCFactorSetReuseFill()`

675:   Level: intermediate

677: .seealso: [](ch_ksp), `PCLU`, `PCCHOLESKY`, `PCILU`, `PCICC`, `PCFactorSetReuseOrdering()`, `PCFactorSetFill()`
678: @*/
679: PetscErrorCode PCFactorSetReuseFill(PC pc, PetscBool flag)
680: {
681:   PetscFunctionBegin;
684:   PetscTryMethod(pc, "PCFactorSetReuseFill_C", (PC, PetscBool), (pc, flag));
685:   PetscFunctionReturn(PETSC_SUCCESS);
686: }

688: PetscErrorCode PCFactorInitialize(PC pc, MatFactorType ftype)
689: {
690:   PC_Factor *fact = (PC_Factor *)pc->data;

692:   PetscFunctionBegin;
693:   PetscCall(MatFactorInfoInitialize(&fact->info));
694:   fact->factortype           = ftype;
695:   fact->info.shifttype       = (PetscReal)MAT_SHIFT_NONE;
696:   fact->info.shiftamount     = 100.0 * PETSC_MACHINE_EPSILON;
697:   fact->info.zeropivot       = 100.0 * PETSC_MACHINE_EPSILON;
698:   fact->info.pivotinblocks   = 1.0;
699:   pc->ops->getfactoredmatrix = PCFactorGetMatrix_Factor;

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

724: PetscErrorCode PCFactorClearComposedFunctions(PC pc)
725: {
726:   PetscFunctionBegin;
727:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetZeroPivot_C", NULL));
728:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorGetZeroPivot_C", NULL));
729:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetShiftType_C", NULL));
730:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorGetShiftType_C", NULL));
731:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetShiftAmount_C", NULL));
732:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorGetShiftAmount_C", NULL));
733:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorGetMatSolverType_C", NULL));
734:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetMatSolverType_C", NULL));
735:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetUpMatSolverType_C", NULL));
736:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetFill_C", NULL));
737:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetMatOrderingType_C", NULL));
738:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetLevels_C", NULL));
739:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorGetLevels_C", NULL));
740:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetAllowDiagonalFill_C", NULL));
741:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorGetAllowDiagonalFill_C", NULL));
742:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetPivotInBlocks_C", NULL));
743:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetUseInPlace_C", NULL));
744:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorGetUseInPlace_C", NULL));
745:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetReuseOrdering_C", NULL));
746:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetReuseFill_C", NULL));
747:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorReorderForNonzeroDiagonal_C", NULL));
748:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetDropTolerance_C", NULL));
749:   PetscFunctionReturn(PETSC_SUCCESS);
750: }