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;

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

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

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

 45:   lu->reuseordering = flag;
 46:   return 0;
 47: }

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

 53:   lu->reusefill = flag;
 54:   return 0;
 55: }

 57: static PetscErrorCode  PCFactorSetUseInPlace_Factor(PC pc,PetscBool flg)
 58: {
 59:   PC_Factor *dir = (PC_Factor*)pc->data;

 61:   dir->inplace = flg;
 62:   return 0;
 63: }

 65: static PetscErrorCode  PCFactorGetUseInPlace_Factor(PC pc,PetscBool *flg)
 66: {
 67:   PC_Factor *dir = (PC_Factor*)pc->data;

 69:   *flg = dir->inplace;
 70:   return 0;
 71: }

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

 77:   Input Parameter:
 78: .  pc  - the preconditioner context

 80:   Notes:
 81:     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.

 83:   Level: intermediate

 85: .seealso: PCFactorSetMatSolverType(), PCFactorGetMatrix()
 86: @*/
 87: PetscErrorCode PCFactorSetUpMatSolverType(PC pc)
 88: {
 90:   PetscTryMethod(pc,"PCFactorSetUpMatSolverType_C",(PC),(pc));
 91:   return 0;
 92: }

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

 97:    Logically Collective on PC

 99:    Input Parameters:
100: +  pc - the preconditioner context
101: -  zero - all pivots smaller than this will be considered zero

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

106:    Level: intermediate

108: .seealso: PCFactorSetShiftType(), PCFactorSetShiftAmount()
109: @*/
110: PetscErrorCode  PCFactorSetZeroPivot(PC pc,PetscReal zero)
111: {
114:   PetscTryMethod(pc,"PCFactorSetZeroPivot_C",(PC,PetscReal),(pc,zero));
115:   return 0;
116: }

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

122:    Logically Collective on PC

124:    Input Parameters:
125: +  pc - the preconditioner context
126: -  shifttype - type of shift; one of MAT_SHIFT_NONE, MAT_SHIFT_NONZERO,  MAT_SHIFT_POSITIVE_DEFINITE, MAT_SHIFT_INBLOCKS

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

131:    Level: intermediate

133: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftAmount()
134: @*/
135: PetscErrorCode  PCFactorSetShiftType(PC pc,MatFactorShiftType shifttype)
136: {
139:   PetscTryMethod(pc,"PCFactorSetShiftType_C",(PC,MatFactorShiftType),(pc,shifttype));
140:   return 0;
141: }

143: /*@
144:    PCFactorSetShiftAmount - adds a quantity to the diagonal of the matrix during
145:      numerical factorization, thus the matrix has nonzero pivots

147:    Logically Collective on PC

149:    Input Parameters:
150: +  pc - the preconditioner context
151: -  shiftamount - amount of shift

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

156:    Level: intermediate

158: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftType()
159: @*/
160: PetscErrorCode  PCFactorSetShiftAmount(PC pc,PetscReal shiftamount)
161: {
164:   PetscTryMethod(pc,"PCFactorSetShiftAmount_C",(PC,PetscReal),(pc,shiftamount));
165:   return 0;
166: }

168: /*@
169:    PCFactorSetDropTolerance - The preconditioner will use an ILU
170:    based on a drop tolerance. (Under development)

172:    Logically Collective on PC

174:    Input Parameters:
175: +  pc - the preconditioner context
176: .  dt - the drop tolerance, try from 1.e-10 to .1
177: .  dtcol - tolerance for column pivot, good values [0.1 to 0.01]
178: -  maxrowcount - the max number of nonzeros allowed in a row, best value
179:                  depends on the number of nonzeros in row of original matrix

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

184:    Level: intermediate

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

189: @*/
190: PetscErrorCode  PCFactorSetDropTolerance(PC pc,PetscReal dt,PetscReal dtcol,PetscInt maxrowcount)
191: {
195:   PetscTryMethod(pc,"PCFactorSetDropTolerance_C",(PC,PetscReal,PetscReal,PetscInt),(pc,dt,dtcol,maxrowcount));
196:   return 0;
197: }

199: /*@
200:    PCFactorGetZeroPivot - Gets the tolerance used to define a zero privot

202:    Not Collective

204:    Input Parameters:
205: .  pc - the preconditioner context

207:    Output Parameter:
208: .  pivot - the tolerance

210:    Level: intermediate

212: .seealso: PCFactorSetZeroPivot()
213: @*/
214: PetscErrorCode  PCFactorGetZeroPivot(PC pc,PetscReal *pivot)
215: {
217:   PetscUseMethod(pc,"PCFactorGetZeroPivot_C",(PC,PetscReal*),(pc,pivot));
218:   return 0;
219: }

221: /*@
222:    PCFactorGetShiftAmount - Gets the tolerance used to define a zero privot

224:    Not Collective

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

229:    Output Parameter:
230: .  shift - how much to shift the diagonal entry

232:    Level: intermediate

234: .seealso: PCFactorSetShiftAmount(), PCFactorSetShiftType(), PCFactorGetShiftType()
235: @*/
236: PetscErrorCode  PCFactorGetShiftAmount(PC pc,PetscReal *shift)
237: {
239:   PetscUseMethod(pc,"PCFactorGetShiftAmount_C",(PC,PetscReal*),(pc,shift));
240:   return 0;
241: }

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

246:    Not Collective

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

251:    Output Parameter:
252: .  type - one of MAT_SHIFT_NONE, MAT_SHIFT_NONZERO,  MAT_SHIFT_POSITIVE_DEFINITE, or MAT_SHIFT_INBLOCKS

254:    Level: intermediate

256: .seealso: PCFactorSetShiftType(), MatFactorShiftType, PCFactorSetShiftAmount(), PCFactorGetShiftAmount()
257: @*/
258: PetscErrorCode  PCFactorGetShiftType(PC pc,MatFactorShiftType *type)
259: {
261:   PetscUseMethod(pc,"PCFactorGetShiftType_C",(PC,MatFactorShiftType*),(pc,type));
262:   return 0;
263: }

265: /*@
266:    PCFactorGetLevels - Gets the number of levels of fill to use.

268:    Logically Collective on PC

270:    Input Parameters:
271: .  pc - the preconditioner context

273:    Output Parameter:
274: .  levels - number of levels of fill

276:    Level: intermediate

278: @*/
279: PetscErrorCode  PCFactorGetLevels(PC pc,PetscInt *levels)
280: {
282:   PetscUseMethod(pc,"PCFactorGetLevels_C",(PC,PetscInt*),(pc,levels));
283:   return 0;
284: }

286: /*@
287:    PCFactorSetLevels - Sets the number of levels of fill to use.

289:    Logically Collective on PC

291:    Input Parameters:
292: +  pc - the preconditioner context
293: -  levels - number of levels of fill

295:    Options Database Key:
296: .  -pc_factor_levels <levels> - Sets fill level

298:    Level: intermediate

300: @*/
301: PetscErrorCode  PCFactorSetLevels(PC pc,PetscInt levels)
302: {
306:   PetscTryMethod(pc,"PCFactorSetLevels_C",(PC,PetscInt),(pc,levels));
307:   return 0;
308: }

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

314:    Logically Collective on PC

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

320:    Options Database Key:
321: .  -pc_factor_diagonal_fill <bool> - allow the diagonal fill

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

326:    Level: intermediate

328: .seealso: PCFactorGetAllowDiagonalFill()
329: @*/
330: PetscErrorCode  PCFactorSetAllowDiagonalFill(PC pc,PetscBool flg)
331: {
333:   PetscTryMethod(pc,"PCFactorSetAllowDiagonalFill_C",(PC,PetscBool),(pc,flg));
334:   return 0;
335: }

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

341:    Logically Collective on PC

343:    Input Parameter:
344: .  pc - the preconditioner context

346:    Output Parameter:
347: .   flg - PETSC_TRUE to turn on, PETSC_FALSE to turn off

349:    Notes:
350:    Does not apply with 0 fill.

352:    Level: intermediate

354: .seealso: PCFactorSetAllowDiagonalFill()
355: @*/
356: PetscErrorCode  PCFactorGetAllowDiagonalFill(PC pc,PetscBool *flg)
357: {
359:   PetscUseMethod(pc,"PCFactorGetAllowDiagonalFill_C",(PC,PetscBool*),(pc,flg));
360:   return 0;
361: }

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

366:    Logically Collective on PC

368:    Input Parameters:
369: +  pc - the preconditioner context
370: -  tol - diagonal entries smaller than this in absolute value are considered zero

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

375:    Level: intermediate

377: .seealso: PCFactorSetFill(), PCFactorSetShiftNonzero(), PCFactorSetZeroPivot(), MatReorderForNonzeroDiagonal()
378: @*/
379: PetscErrorCode  PCFactorReorderForNonzeroDiagonal(PC pc,PetscReal rtol)
380: {
383:   PetscTryMethod(pc,"PCFactorReorderForNonzeroDiagonal_C",(PC,PetscReal),(pc,rtol));
384:   return 0;
385: }

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

390:    Logically Collective on PC

392:    Input Parameters:
393: +  pc - the preconditioner context
394: -  stype - for example, superlu, superlu_dist

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

399:    Level: intermediate

401:    Note:
402:      By default this will use the PETSc factorization if it exists

404: .seealso: MatGetFactor(), MatSolverType, PCFactorGetMatSolverType()
405: @*/
406: PetscErrorCode  PCFactorSetMatSolverType(PC pc,MatSolverType stype)
407: {
409:   PetscTryMethod(pc,"PCFactorSetMatSolverType_C",(PC,MatSolverType),(pc,stype));
410:   return 0;
411: }

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

416:    Not Collective

418:    Input Parameter:
419: .  pc - the preconditioner context

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

424:    Level: intermediate

426: .seealso: MatGetFactor(), MatSolverType, PCFactorGetMatSolverType()
427: @*/
428: PetscErrorCode  PCFactorGetMatSolverType(PC pc,MatSolverType *stype)
429: {
430:   PetscErrorCode (*f)(PC,MatSolverType*);

434:   PetscObjectQueryFunction((PetscObject)pc,"PCFactorGetMatSolverType_C",&f);
435:   if (f) (*f)(pc,stype);
436:   else *stype = NULL;
437:   return 0;
438: }

440: /*@
441:    PCFactorSetFill - Indicate the amount of fill you expect in the factored matrix,
442:    fill = number nonzeros in factor/number nonzeros in original matrix.

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

446:    Input Parameters:
447: +  pc - the preconditioner context
448: -  fill - amount of expected fill

450:    Options Database Key:
451: .  -pc_factor_fill <fill> - Sets fill amount

453:    Level: intermediate

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

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

463: @*/
464: PetscErrorCode  PCFactorSetFill(PC pc,PetscReal fill)
465: {
468:   PetscTryMethod(pc,"PCFactorSetFill_C",(PC,PetscReal),(pc,fill));
469:   return 0;
470: }

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

480:    Logically Collective on PC

482:    Input Parameters:
483: +  pc - the preconditioner context
484: -  flg - PETSC_TRUE to enable, PETSC_FALSE to disable

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

489:    Notes:
490:    PCFactorSetUseInplace() can only be used with the KSP method KSPPREONLY or when
491:    a different matrix is provided for the multiply and the preconditioner in
492:    a call to KSPSetOperators().
493:    This is because the Krylov space methods require an application of the
494:    matrix multiplication, which is not possible here because the matrix has
495:    been factored in-place, replacing the original matrix.

497:    Level: intermediate

499: .seealso: PCFactorGetUseInPlace()
500: @*/
501: PetscErrorCode  PCFactorSetUseInPlace(PC pc,PetscBool flg)
502: {
504:   PetscTryMethod(pc,"PCFactorSetUseInPlace_C",(PC,PetscBool),(pc,flg));
505:   return 0;
506: }

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

511:    Logically Collective on PC

513:    Input Parameter:
514: .  pc - the preconditioner context

516:    Output Parameter:
517: .  flg - PETSC_TRUE to enable, PETSC_FALSE to disable

519:    Level: intermediate

521: .seealso: PCFactorSetUseInPlace()
522: @*/
523: PetscErrorCode  PCFactorGetUseInPlace(PC pc,PetscBool *flg)
524: {
526:   PetscUseMethod(pc,"PCFactorGetUseInPlace_C",(PC,PetscBool*),(pc,flg));
527:   return 0;
528: }

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

534:     Logically Collective on PC

536:     Input Parameters:
537: +   pc - the preconditioner context
538: -   ordering - the matrix ordering name, for example, MATORDERINGND or MATORDERINGRCM

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

543:     Level: intermediate

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

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

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

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

555: .seealso: MatOrderingType

557: @*/
558: PetscErrorCode  PCFactorSetMatOrderingType(PC pc,MatOrderingType ordering)
559: {
561:   PetscTryMethod(pc,"PCFactorSetMatOrderingType_C",(PC,MatOrderingType),(pc,ordering));
562:   return 0;
563: }

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

570:     Logically Collective on PC

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

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

579:     Level: intermediate

581: .seealso: PCILUSetMatOrdering(), PCFactorSetPivotInBlocks()
582: @*/
583: PetscErrorCode  PCFactorSetColumnPivot(PC pc,PetscReal dtcol)
584: {
587:   PetscTryMethod(pc,"PCFactorSetColumnPivot_C",(PC,PetscReal),(pc,dtcol));
588:   return 0;
589: }

591: /*@
592:     PCFactorSetPivotInBlocks - Determines if pivoting is done while factoring each block
593:       with BAIJ or SBAIJ matrices

595:     Logically Collective on PC

597:     Input Parameters:
598: +   pc - the preconditioner context
599: -   pivot - PETSC_TRUE or PETSC_FALSE

601:     Options Database Key:
602: .   -pc_factor_pivot_in_blocks <true,false> - Pivot inside matrix dense blocks for BAIJ and SBAIJ

604:     Level: intermediate

606: .seealso: PCILUSetMatOrdering(), PCFactorSetColumnPivot()
607: @*/
608: PetscErrorCode  PCFactorSetPivotInBlocks(PC pc,PetscBool pivot)
609: {
612:   PetscTryMethod(pc,"PCFactorSetPivotInBlocks_C",(PC,PetscBool),(pc,pivot));
613:   return 0;
614: }

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

620:    Logically Collective on PC

622:    Input Parameters:
623: +  pc - the preconditioner context
624: -  flag - PETSC_TRUE to reuse else PETSC_FALSE

626:    Options Database Key:
627: .  -pc_factor_reuse_fill - Activates PCFactorSetReuseFill()

629:    Level: intermediate

631: .seealso: PCFactorSetReuseOrdering()
632: @*/
633: PetscErrorCode  PCFactorSetReuseFill(PC pc,PetscBool flag)
634: {
637:   PetscTryMethod(pc,"PCFactorSetReuseFill_C",(PC,PetscBool),(pc,flag));
638:   return 0;
639: }

641: PetscErrorCode PCFactorInitialize(PC pc,MatFactorType ftype)
642: {
643:   PC_Factor      *fact = (PC_Factor*)pc->data;

645:   MatFactorInfoInitialize(&fact->info);
646:   fact->factortype           = ftype;
647:   fact->info.shifttype       = (PetscReal)MAT_SHIFT_NONE;
648:   fact->info.shiftamount     = 100.0*PETSC_MACHINE_EPSILON;
649:   fact->info.zeropivot       = 100.0*PETSC_MACHINE_EPSILON;
650:   fact->info.pivotinblocks   = 1.0;
651:   pc->ops->getfactoredmatrix = PCFactorGetMatrix_Factor;

653:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetZeroPivot_C",PCFactorSetZeroPivot_Factor);
654:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetZeroPivot_C",PCFactorGetZeroPivot_Factor);
655:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftType_C",PCFactorSetShiftType_Factor);
656:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetShiftType_C",PCFactorGetShiftType_Factor);
657:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftAmount_C",PCFactorSetShiftAmount_Factor);
658:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetShiftAmount_C",PCFactorGetShiftAmount_Factor);
659:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetMatSolverType_C",PCFactorGetMatSolverType_Factor);
660:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatSolverType_C",PCFactorSetMatSolverType_Factor);
661:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUpMatSolverType_C",PCFactorSetUpMatSolverType_Factor);
662:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetFill_C",PCFactorSetFill_Factor);
663:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatOrderingType_C",PCFactorSetMatOrderingType_Factor);
664:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetLevels_C",PCFactorSetLevels_Factor);
665:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetLevels_C",PCFactorGetLevels_Factor);
666:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetAllowDiagonalFill_C",PCFactorSetAllowDiagonalFill_Factor);
667:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetAllowDiagonalFill_C",PCFactorGetAllowDiagonalFill_Factor);
668:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetPivotInBlocks_C",PCFactorSetPivotInBlocks_Factor);
669:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUseInPlace_C",PCFactorSetUseInPlace_Factor);
670:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetUseInPlace_C",PCFactorGetUseInPlace_Factor);
671:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseOrdering_C",PCFactorSetReuseOrdering_Factor);
672:   PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseFill_C",PCFactorSetReuseFill_Factor);
673:   return 0;
674: }