Actual source code: factor.c


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

  5: /*
  6:     If an ordering is not yet set and the matrix is available determine a default ordering
  7: */
  8: PetscErrorCode PCFactorSetDefaultOrdering_Factor(PC pc)
  9: {
 10:   Mat            B;
 11:   PetscBool      foundmtype,flg;

 15:   if (pc->pmat) {
 16:     PC_Factor *fact = (PC_Factor*)pc->data;
 17:     MatSolverTypeGet(fact->solvertype,((PetscObject)pc->pmat)->type_name,fact->factortype,NULL,&foundmtype,NULL);
 18:     if (foundmtype) {
 19:       if (!fact->fact) {
 20:         MatGetFactor(pc->pmat,fact->solvertype,fact->factortype,&fact->fact);
 21:       } else if (!fact->fact->assembled) {
 22:         PetscStrcmp(fact->solvertype,fact->fact->solvertype,&flg);
 23:         if (!flg) {
 24:           MatGetFactor(pc->pmat,fact->solvertype,fact->factortype,&B);
 25:           MatHeaderReplace(fact->fact,&B);
 26:         }
 27:       }
 28:       if (!fact->ordering) {
 29:         PetscBool       canuseordering;
 30:         MatOrderingType otype;

 32:         MatFactorGetCanUseOrdering(fact->fact,&canuseordering);
 33:         if (canuseordering) {
 34:           MatFactorGetPreferredOrdering(fact->fact,fact->factortype,&otype);
 35:         } else otype = MATORDERINGEXTERNAL;
 36:         PetscStrallocpy(otype,(char **)&fact->ordering);
 37:       }
 38:     }
 39:   }
 40:   return(0);
 41: }

 43: static PetscErrorCode PCFactorSetReuseOrdering_Factor(PC pc,PetscBool flag)
 44: {
 45:   PC_Factor *lu = (PC_Factor*)pc->data;

 48:   lu->reuseordering = flag;
 49:   return(0);
 50: }

 52: static PetscErrorCode PCFactorSetReuseFill_Factor(PC pc,PetscBool flag)
 53: {
 54:   PC_Factor *lu = (PC_Factor*)pc->data;

 57:   lu->reusefill = flag;
 58:   return(0);
 59: }

 61: static PetscErrorCode  PCFactorSetUseInPlace_Factor(PC pc,PetscBool flg)
 62: {
 63:   PC_Factor *dir = (PC_Factor*)pc->data;

 66:   dir->inplace = flg;
 67:   return(0);
 68: }

 70: static PetscErrorCode  PCFactorGetUseInPlace_Factor(PC pc,PetscBool *flg)
 71: {
 72:   PC_Factor *dir = (PC_Factor*)pc->data;

 75:   *flg = dir->inplace;
 76:   return(0);
 77: }

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

 83:   Input Parameter:
 84: .  pc  - the preconditioner context

 86:   Notes:
 87:     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.

 89:   Level: intermediate

 91: .seealso: PCFactorSetMatSolverType(), PCFactorGetMatrix()
 92: @*/
 93: PetscErrorCode PCFactorSetUpMatSolverType(PC pc)
 94: {

 99:   PetscTryMethod(pc,"PCFactorSetUpMatSolverType_C",(PC),(pc));
100:   return(0);
101: }

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

106:    Logically Collective on PC

108:    Input Parameters:
109: +  pc - the preconditioner context
110: -  zero - all pivots smaller than this will be considered zero

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

115:    Level: intermediate

117: .seealso: PCFactorSetShiftType(), PCFactorSetShiftAmount()
118: @*/
119: PetscErrorCode  PCFactorSetZeroPivot(PC pc,PetscReal zero)
120: {

126:   PetscTryMethod(pc,"PCFactorSetZeroPivot_C",(PC,PetscReal),(pc,zero));
127:   return(0);
128: }

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

134:    Logically Collective on PC

136:    Input Parameters:
137: +  pc - the preconditioner context
138: -  shifttype - type of shift; one of MAT_SHIFT_NONE, MAT_SHIFT_NONZERO,  MAT_SHIFT_POSITIVE_DEFINITE, MAT_SHIFT_INBLOCKS

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

143:    Level: intermediate

145: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftAmount()
146: @*/
147: PetscErrorCode  PCFactorSetShiftType(PC pc,MatFactorShiftType shifttype)
148: {

154:   PetscTryMethod(pc,"PCFactorSetShiftType_C",(PC,MatFactorShiftType),(pc,shifttype));
155:   return(0);
156: }

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

162:    Logically Collective on PC

164:    Input Parameters:
165: +  pc - the preconditioner context
166: -  shiftamount - amount of shift

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

171:    Level: intermediate

173: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftType()
174: @*/
175: PetscErrorCode  PCFactorSetShiftAmount(PC pc,PetscReal shiftamount)
176: {

182:   PetscTryMethod(pc,"PCFactorSetShiftAmount_C",(PC,PetscReal),(pc,shiftamount));
183:   return(0);
184: }

186: /*@
187:    PCFactorSetDropTolerance - The preconditioner will use an ILU
188:    based on a drop tolerance. (Under development)

190:    Logically Collective on PC

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:       There are NO default values for the 3 parameters, you must set them with reasonable values for your
205:       matrix. We don't know how to compute reasonable values.

207: @*/
208: PetscErrorCode  PCFactorSetDropTolerance(PC pc,PetscReal dt,PetscReal dtcol,PetscInt maxrowcount)
209: {

216:   PetscTryMethod(pc,"PCFactorSetDropTolerance_C",(PC,PetscReal,PetscReal,PetscInt),(pc,dt,dtcol,maxrowcount));
217:   return(0);
218: }

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

223:    Not Collective

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

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

231:    Level: intermediate

233: .seealso: PCFactorSetZeroPivot()
234: @*/
235: PetscErrorCode  PCFactorGetZeroPivot(PC pc,PetscReal *pivot)
236: {

241:   PetscUseMethod(pc,"PCFactorGetZeroPivot_C",(PC,PetscReal*),(pc,pivot));
242:   return(0);
243: }

245: /*@
246:    PCFactorGetShiftAmount - Gets the tolerance used to define a zero privot

248:    Not Collective

250:    Input Parameters:
251: .  pc - the preconditioner context

253:    Output Parameter:
254: .  shift - how much to shift the diagonal entry

256:    Level: intermediate

258: .seealso: PCFactorSetShiftAmount(), PCFactorSetShiftType(), PCFactorGetShiftType()
259: @*/
260: PetscErrorCode  PCFactorGetShiftAmount(PC pc,PetscReal *shift)
261: {

266:   PetscUseMethod(pc,"PCFactorGetShiftAmount_C",(PC,PetscReal*),(pc,shift));
267:   return(0);
268: }

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

273:    Not Collective

275:    Input Parameters:
276: .  pc - the preconditioner context

278:    Output Parameter:
279: .  type - one of MAT_SHIFT_NONE, MAT_SHIFT_NONZERO,  MAT_SHIFT_POSITIVE_DEFINITE, or MAT_SHIFT_INBLOCKS

281:    Level: intermediate

283: .seealso: PCFactorSetShiftType(), MatFactorShiftType, PCFactorSetShiftAmount(), PCFactorGetShiftAmount()
284: @*/
285: PetscErrorCode  PCFactorGetShiftType(PC pc,MatFactorShiftType *type)
286: {

291:   PetscUseMethod(pc,"PCFactorGetShiftType_C",(PC,MatFactorShiftType*),(pc,type));
292:   return(0);
293: }

295: /*@
296:    PCFactorGetLevels - Gets the number of levels of fill to use.

298:    Logically Collective on PC

300:    Input Parameters:
301: .  pc - the preconditioner context

303:    Output Parameter:
304: .  levels - number of levels of fill

306:    Level: intermediate

308: @*/
309: PetscErrorCode  PCFactorGetLevels(PC pc,PetscInt *levels)
310: {

315:   PetscUseMethod(pc,"PCFactorGetLevels_C",(PC,PetscInt*),(pc,levels));
316:   return(0);
317: }

319: /*@
320:    PCFactorSetLevels - Sets the number of levels of fill to use.

322:    Logically Collective on PC

324:    Input Parameters:
325: +  pc - the preconditioner context
326: -  levels - number of levels of fill

328:    Options Database Key:
329: .  -pc_factor_levels <levels> - Sets fill level

331:    Level: intermediate

333: @*/
334: PetscErrorCode  PCFactorSetLevels(PC pc,PetscInt levels)
335: {

340:   if (levels < 0) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"negative levels");
342:   PetscTryMethod(pc,"PCFactorSetLevels_C",(PC,PetscInt),(pc,levels));
343:   return(0);
344: }

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

350:    Logically Collective on PC

352:    Input Parameters:
353: +  pc - the preconditioner context
354: -  flg - PETSC_TRUE to turn on, PETSC_FALSE to turn off

356:    Options Database Key:
357: .  -pc_factor_diagonal_fill

359:    Notes:
360:    Does not apply with 0 fill.

362:    Level: intermediate

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

372:   PetscTryMethod(pc,"PCFactorSetAllowDiagonalFill_C",(PC,PetscBool),(pc,flg));
373:   return(0);
374: }

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

380:    Logically Collective on PC

382:    Input Parameter:
383: .  pc - the preconditioner context

385:    Output Parameter:
386: .   flg - PETSC_TRUE to turn on, PETSC_FALSE to turn off

388:    Options Database Key:
389: .  -pc_factor_diagonal_fill

391:    Notes:
392:    Does not apply with 0 fill.

394:    Level: intermediate

396: .seealso: PCFactorSetAllowDiagonalFill()
397: @*/
398: PetscErrorCode  PCFactorGetAllowDiagonalFill(PC pc,PetscBool *flg)
399: {

404:   PetscUseMethod(pc,"PCFactorGetAllowDiagonalFill_C",(PC,PetscBool*),(pc,flg));
405:   return(0);
406: }

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

411:    Logically Collective on PC

413:    Input Parameters:
414: +  pc - the preconditioner context
415: -  tol - diagonal entries smaller than this in absolute value are considered zero

417:    Options Database Key:
418: .  -pc_factor_nonzeros_along_diagonal <tol>

420:    Level: intermediate

422: .seealso: PCFactorSetFill(), PCFactorSetShiftNonzero(), PCFactorSetZeroPivot(), MatReorderForNonzeroDiagonal()
423: @*/
424: PetscErrorCode  PCFactorReorderForNonzeroDiagonal(PC pc,PetscReal rtol)
425: {

431:   PetscTryMethod(pc,"PCFactorReorderForNonzeroDiagonal_C",(PC,PetscReal),(pc,rtol));
432:   return(0);
433: }

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

438:    Logically Collective on PC

440:    Input Parameters:
441: +  pc - the preconditioner context
442: -  stype - for example, superlu, superlu_dist

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

447:    Level: intermediate

449:    Note:
450:      By default this will use the PETSc factorization if it exists

452: .seealso: MatGetFactor(), MatSolverType, PCFactorGetMatSolverType()
453: @*/
454: PetscErrorCode  PCFactorSetMatSolverType(PC pc,MatSolverType stype)
455: {

460:   PetscTryMethod(pc,"PCFactorSetMatSolverType_C",(PC,MatSolverType),(pc,stype));
461:   return(0);
462: }

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

467:    Not Collective

469:    Input Parameter:
470: .  pc - the preconditioner context

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

475:    Level: intermediate

477: .seealso: MatGetFactor(), MatSolverType, PCFactorGetMatSolverType()
478: @*/
479: PetscErrorCode  PCFactorGetMatSolverType(PC pc,MatSolverType *stype)
480: {
481:   PetscErrorCode ierr,(*f)(PC,MatSolverType*);

485:   PetscObjectQueryFunction((PetscObject)pc,"PCFactorGetMatSolverType_C",&f);
486:   if (f) {
487:     (*f)(pc,stype);
488:   } else {
489:     *stype = NULL;
490:   }
491:   return(0);
492: }

494: /*@
495:    PCFactorSetFill - Indicate the amount of fill you expect in the factored matrix,
496:    fill = number nonzeros in factor/number nonzeros in original matrix.

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

500:    Input Parameters:
501: +  pc - the preconditioner context
502: -  fill - amount of expected fill

504:    Options Database Key:
505: .  -pc_factor_fill <fill> - Sets fill amount

507:    Level: intermediate

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

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

517: @*/
518: PetscErrorCode  PCFactorSetFill(PC pc,PetscReal fill)
519: {

524:   if (fill < 1.0) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Fill factor cannot be less then 1.0");
525:   PetscTryMethod(pc,"PCFactorSetFill_C",(PC,PetscReal),(pc,fill));
526:   return(0);
527: }

529: /*@
530:    PCFactorSetUseInPlace - Tells the system to do an in-place factorization.
531:    For dense matrices, this enables the solution of much larger problems.
532:    For sparse matrices the factorization cannot be done truly in-place
533:    so this does not save memory during the factorization, but after the matrix
534:    is factored, the original unfactored matrix is freed, thus recovering that
535:    space. For ICC(0) and ILU(0) with the default natural ordering the factorization is done efficiently in-place.

537:    Logically Collective on PC

539:    Input Parameters:
540: +  pc - the preconditioner context
541: -  flg - PETSC_TRUE to enable, PETSC_FALSE to disable

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

546:    Notes:
547:    PCFactorSetUseInplace() can only be used with the KSP method KSPPREONLY or when
548:    a different matrix is provided for the multiply and the preconditioner in
549:    a call to KSPSetOperators().
550:    This is because the Krylov space methods require an application of the
551:    matrix multiplication, which is not possible here because the matrix has
552:    been factored in-place, replacing the original matrix.

554:    Level: intermediate

556: .seealso: PCFactorGetUseInPlace()
557: @*/
558: PetscErrorCode  PCFactorSetUseInPlace(PC pc,PetscBool flg)
559: {

564:   PetscTryMethod(pc,"PCFactorSetUseInPlace_C",(PC,PetscBool),(pc,flg));
565:   return(0);
566: }

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

571:    Logically Collective on PC

573:    Input Parameter:
574: .  pc - the preconditioner context

576:    Output Parameter:
577: .  flg - PETSC_TRUE to enable, PETSC_FALSE to disable

579:    Level: intermediate

581: .seealso: PCFactorSetUseInPlace()
582: @*/
583: PetscErrorCode  PCFactorGetUseInPlace(PC pc,PetscBool *flg)
584: {

589:   PetscUseMethod(pc,"PCFactorGetUseInPlace_C",(PC,PetscBool*),(pc,flg));
590:   return(0);
591: }

593: /*@C
594:     PCFactorSetMatOrderingType - Sets the ordering routine (to reduce fill) to
595:     be used in the LU, ILU, Cholesky, and ICC factorizations.

597:     Logically Collective on PC

599:     Input Parameters:
600: +   pc - the preconditioner context
601: -   ordering - the matrix ordering name, for example, MATORDERINGND or MATORDERINGRCM

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

606:     Level: intermediate

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

611:      For Cholesky and ICC and the SBAIJ format the only reordering available is natural since only the upper half of the matrix is stored
612:      and reordering this matrix is very expensive.

614:       You can use a SeqAIJ matrix with Cholesky and ICC and use any ordering.

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

618: .seealso: MatOrderingType

620: @*/
621: PetscErrorCode  PCFactorSetMatOrderingType(PC pc,MatOrderingType ordering)
622: {

627:   PetscTryMethod(pc,"PCFactorSetMatOrderingType_C",(PC,MatOrderingType),(pc,ordering));
628:   return(0);
629: }

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

636:     Logically Collective on PC

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

642:     Options Database Key:
643: .   -pc_factor_pivoting <dtcol>

645:     Level: intermediate

647: .seealso: PCILUSetMatOrdering(), PCFactorSetPivotInBlocks()
648: @*/
649: PetscErrorCode  PCFactorSetColumnPivot(PC pc,PetscReal dtcol)
650: {

656:   PetscTryMethod(pc,"PCFactorSetColumnPivot_C",(PC,PetscReal),(pc,dtcol));
657:   return(0);
658: }

660: /*@
661:     PCFactorSetPivotInBlocks - Determines if pivoting is done while factoring each block
662:       with BAIJ or SBAIJ matrices

664:     Logically Collective on PC

666:     Input Parameters:
667: +   pc - the preconditioner context
668: -   pivot - PETSC_TRUE or PETSC_FALSE

670:     Options Database Key:
671: .   -pc_factor_pivot_in_blocks <true,false>

673:     Level: intermediate

675: .seealso: PCILUSetMatOrdering(), PCFactorSetColumnPivot()
676: @*/
677: PetscErrorCode  PCFactorSetPivotInBlocks(PC pc,PetscBool pivot)
678: {

684:   PetscTryMethod(pc,"PCFactorSetPivotInBlocks_C",(PC,PetscBool),(pc,pivot));
685:   return(0);
686: }

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

692:    Logically Collective on PC

694:    Input Parameters:
695: +  pc - the preconditioner context
696: -  flag - PETSC_TRUE to reuse else PETSC_FALSE

698:    Options Database Key:
699: .  -pc_factor_reuse_fill - Activates PCFactorSetReuseFill()

701:    Level: intermediate

703: .seealso: PCFactorSetReuseOrdering()
704: @*/
705: PetscErrorCode  PCFactorSetReuseFill(PC pc,PetscBool flag)
706: {

712:   PetscTryMethod(pc,"PCFactorSetReuseFill_C",(PC,PetscBool),(pc,flag));
713:   return(0);
714: }

716: PetscErrorCode PCFactorInitialize(PC pc,MatFactorType ftype)
717: {
719:   PC_Factor      *fact = (PC_Factor*)pc->data;

722:   MatFactorInfoInitialize(&fact->info);
723:   fact->factortype           = ftype;
724:   fact->info.shifttype       = (PetscReal)MAT_SHIFT_NONE;
725:   fact->info.shiftamount     = 100.0*PETSC_MACHINE_EPSILON;
726:   fact->info.zeropivot       = 100.0*PETSC_MACHINE_EPSILON;
727:   fact->info.pivotinblocks   = 1.0;
728:   pc->ops->getfactoredmatrix = PCFactorGetMatrix_Factor;

730:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetZeroPivot_C",PCFactorSetZeroPivot_Factor);
731:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetZeroPivot_C",PCFactorGetZeroPivot_Factor);
732:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftType_C",PCFactorSetShiftType_Factor);
733:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetShiftType_C",PCFactorGetShiftType_Factor);
734:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftAmount_C",PCFactorSetShiftAmount_Factor);
735:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetShiftAmount_C",PCFactorGetShiftAmount_Factor);
736:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetMatSolverType_C",PCFactorGetMatSolverType_Factor);
737:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatSolverType_C",PCFactorSetMatSolverType_Factor);
738:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUpMatSolverType_C",PCFactorSetUpMatSolverType_Factor);
739:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetFill_C",PCFactorSetFill_Factor);
740:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatOrderingType_C",PCFactorSetMatOrderingType_Factor);
741:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetLevels_C",PCFactorSetLevels_Factor);
742:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetLevels_C",PCFactorGetLevels_Factor);
743:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetAllowDiagonalFill_C",PCFactorSetAllowDiagonalFill_Factor);
744:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetAllowDiagonalFill_C",PCFactorGetAllowDiagonalFill_Factor);
745:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetPivotInBlocks_C",PCFactorSetPivotInBlocks_Factor);
746:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUseInPlace_C",PCFactorSetUseInPlace_Factor);
747:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetUseInPlace_C",PCFactorGetUseInPlace_Factor);
748:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseOrdering_C",PCFactorSetReuseOrdering_Factor);
749:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseFill_C",PCFactorSetReuseFill_Factor);
750:   return(0);
751: }