Actual source code: ilu.c

  1: /*$Id: ilu.c,v 1.165 2001/04/13 03:00:55 bsmith Exp $*/
  2: /*
  3:    Defines a ILU factorization preconditioner for any Mat implementation
  4: */
 5:  #include src/sles/pc/pcimpl.h
 6:  #include src/sles/pc/impls/ilu/ilu.h
 7:  #include src/mat/matimpl.h

  9: /* ------------------------------------------------------------------------------------------*/
 10: EXTERN_C_BEGIN
 11: int PCILUSetDamping_ILU(PC pc,PetscReal damping)
 12: {
 13:   PC_ILU *dir;

 16:   dir = (PC_ILU*)pc->data;
 17:   dir->info.damping = damping;
 18:   dir->info.damp    = 1.0;
 19:   return(0);
 20: }
 21: EXTERN_C_END

 23: EXTERN_C_BEGIN
 24: int PCILUSetUseDropTolerance_ILU(PC pc,PetscReal dt,PetscReal dtcol,int dtcount)
 25: {
 26:   PC_ILU *ilu;

 29:   ilu = (PC_ILU*)pc->data;
 30:   ilu->usedt         = PETSC_TRUE;
 31:   ilu->info.dt       = dt;
 32:   ilu->info.dtcol    = dtcol;
 33:   ilu->info.dtcount  = dtcount;
 34:   ilu->info.fill     = PETSC_DEFAULT;
 35:   return(0);
 36: }
 37: EXTERN_C_END

 39: EXTERN_C_BEGIN
 40: int PCILUSetFill_ILU(PC pc,PetscReal fill)
 41: {
 42:   PC_ILU *dir;

 45:   dir            = (PC_ILU*)pc->data;
 46:   dir->info.fill = fill;
 47:   return(0);
 48: }
 49: EXTERN_C_END

 51: EXTERN_C_BEGIN
 52: int PCILUSetMatOrdering_ILU(PC pc,MatOrderingType ordering)
 53: {
 54:   PC_ILU *dir = (PC_ILU*)pc->data;
 55:   int    ierr;
 56: 
 58:   PetscStrfree(dir->ordering);
 59:   PetscStrallocpy(ordering,&dir->ordering);
 60:   return(0);
 61: }
 62: EXTERN_C_END

 64: EXTERN_C_BEGIN
 65: int PCILUSetReuseOrdering_ILU(PC pc,PetscTruth flag)
 66: {
 67:   PC_ILU *ilu;

 70:   ilu                = (PC_ILU*)pc->data;
 71:   ilu->reuseordering = flag;
 72:   return(0);
 73: }
 74: EXTERN_C_END

 76: EXTERN_C_BEGIN
 77: int PCILUDTSetReuseFill_ILUDT(PC pc,PetscTruth flag)
 78: {
 79:   PC_ILU *ilu;

 82:   ilu = (PC_ILU*)pc->data;
 83:   ilu->reusefill = flag;
 84:   if (flag) SETERRQ(1,"Not yet supported");
 85:   return(0);
 86: }
 87: EXTERN_C_END

 89: EXTERN_C_BEGIN
 90: int PCILUSetLevels_ILU(PC pc,int levels)
 91: {
 92:   PC_ILU *ilu;

 95:   ilu = (PC_ILU*)pc->data;
 96:   ilu->info.levels = levels;
 97:   return(0);
 98: }
 99: EXTERN_C_END

101: EXTERN_C_BEGIN
102: int PCILUSetUseInPlace_ILU(PC pc)
103: {
104:   PC_ILU *dir;

107:   dir          = (PC_ILU*)pc->data;
108:   dir->inplace = PETSC_TRUE;
109:   return(0);
110: }
111: EXTERN_C_END

113: EXTERN_C_BEGIN
114: int PCILUSetAllowDiagonalFill_ILU(PC pc)
115: {
116:   PC_ILU *dir;

119:   dir                 = (PC_ILU*)pc->data;
120:   dir->info.diagonal_fill = 1;
121:   return(0);
122: }
123: EXTERN_C_END

125: /* ------------------------------------------------------------------------------------------*/
126: /*@
127:    PCILUSetDamping - adds this quantity to the diagonal of the matrix during the 
128:      ILU numerical factorization

130:    Collective on PC
131:    
132:    Input Parameters:
133: +  pc - the preconditioner context
134: -  damping - amount of damping

136:    Options Database Key:
137: .  -pc_ilu_damping <damping> - Sets damping amount

139:    Level: intermediate

141: .keywords: PC, set, factorization, direct, fill

143: .seealso: PCILUSetFill(), PCLUSetDamp()
144: @*/
145: int PCILUSetDamping(PC pc,PetscReal damping)
146: {
147:   int ierr,(*f)(PC,PetscReal);

151:   PetscObjectQueryFunction((PetscObject)pc,"PCILUSetDamping_C",(void (**)())&f);
152:   if (f) {
153:     (*f)(pc,damping);
154:   }
155:   return(0);
156: }

158: /*@
159:    PCILUSetUseDropTolerance - The preconditioner will use an ILU 
160:    based on a drop tolerance.

162:    Collective on PC

164:    Input Parameters:
165: +  pc - the preconditioner context
166: .  dt - the drop tolerance, try from 1.e-10 to .1
167: .  dtcol - tolerance for column pivot, good values [0.1 to 0.01]
168: -  maxrowcount - the max number of nonzeros allowed in a row, best value
169:                  depends on the number of nonzeros in row of original matrix

171:    Options Database Key:
172: .  -pc_ilu_use_drop_tolerance <dt,dtcol,maxrowcount> - Sets drop tolerance

174:    Level: intermediate

176:     Notes:
177:       This uses the iludt() code of Saad's SPARSKIT package

179: .keywords: PC, levels, reordering, factorization, incomplete, ILU
180: @*/
181: int PCILUSetUseDropTolerance(PC pc,PetscReal dt,PetscReal dtcol,int maxrowcount)
182: {
183:   int ierr,(*f)(PC,PetscReal,PetscReal,int);

187:   PetscObjectQueryFunction((PetscObject)pc,"PCILUSetUseDropTolerance_C",(void (**)())&f);
188:   if (f) {
189:     (*f)(pc,dt,dtcol,maxrowcount);
190:   }
191:   return(0);
192: }

194: /*@
195:    PCILUSetFill - Indicate the amount of fill you expect in the factored matrix,
196:    where fill = number nonzeros in factor/number nonzeros in original matrix.

198:    Collective on PC

200:    Input Parameters:
201: +  pc - the preconditioner context
202: -  fill - amount of expected fill

204:    Options Database Key:
205: $  -pc_ilu_fill <fill>

207:    Note:
208:    For sparse matrix factorizations it is difficult to predict how much 
209:    fill to expect. By running with the option -log_info PETSc will print the 
210:    actual amount of fill used; allowing you to set the value accurately for
211:    future runs. But default PETSc uses a value of 1.0

213:    Level: intermediate

215: .keywords: PC, set, factorization, direct, fill

217: .seealso: PCLUSetFill()
218: @*/
219: int PCILUSetFill(PC pc,PetscReal fill)
220: {
221:   int ierr,(*f)(PC,PetscReal);

225:   if (fill < 1.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Fill factor cannot be less than 1.0");
226:   PetscObjectQueryFunction((PetscObject)pc,"PCILUSetFill_C",(void (**)())&f);
227:   if (f) {
228:     (*f)(pc,fill);
229:   }
230:   return(0);
231: }

233: /*@
234:     PCILUSetMatOrdering - Sets the ordering routine (to reduce fill) to 
235:     be used in the ILU factorization.

237:     Collective on PC

239:     Input Parameters:
240: +   pc - the preconditioner context
241: -   ordering - the matrix ordering name, for example, MATORDERING_ND or MATORDERING_RCM

243:     Options Database Key:
244: .   -pc_ilu_mat_ordering_type <nd,rcm,...> - Sets ordering routine

246:     Level: intermediate

248:     Notes: natural ordering is used by default

250: .seealso: PCLUSetMatOrdering()

252: .keywords: PC, ILU, set, matrix, reordering

254: @*/
255: int PCILUSetMatOrdering(PC pc,MatOrderingType ordering)
256: {
257:   int ierr,(*f)(PC,MatOrderingType);

261:   PetscObjectQueryFunction((PetscObject)pc,"PCILUSetMatOrdering_C",(void (**)())&f);
262:   if (f) {
263:     (*f)(pc,ordering);
264:   }
265:   return(0);
266: }

268: /*@
269:    PCILUSetReuseOrdering - When similar matrices are factored, this
270:    causes the ordering computed in the first factor to be used for all
271:    following factors; applies to both fill and drop tolerance ILUs.

273:    Collective on PC

275:    Input Parameters:
276: +  pc - the preconditioner context
277: -  flag - PETSC_TRUE to reuse else PETSC_FALSE

279:    Options Database Key:
280: .  -pc_ilu_reuse_ordering - Activate PCILUSetReuseOrdering()

282:    Level: intermediate

284: .keywords: PC, levels, reordering, factorization, incomplete, ILU

286: .seealso: PCILUDTSetReuseFill(), PCLUSetReuseOrdering(), PCLUSetReuseFill()
287: @*/
288: int PCILUSetReuseOrdering(PC pc,PetscTruth flag)
289: {
290:   int ierr,(*f)(PC,PetscTruth);

294:   PetscObjectQueryFunction((PetscObject)pc,"PCILUSetReuseOrdering_C",(void (**)())&f);
295:   if (f) {
296:     (*f)(pc,flag);
297:   }
298:   return(0);
299: }

301: /*@
302:    PCILUDTSetReuseFill - When matrices with same nonzero structure are ILUDT factored,
303:    this causes later ones to use the fill computed in the initial factorization.

305:    Collective on PC

307:    Input Parameters:
308: +  pc - the preconditioner context
309: -  flag - PETSC_TRUE to reuse else PETSC_FALSE

311:    Options Database Key:
312: .  -pc_iludt_reuse_fill - Activates PCILUDTSetReuseFill()

314:    Level: intermediate

316: .keywords: PC, levels, reordering, factorization, incomplete, ILU

318: .seealso: PCILUSetReuseOrdering(), PCLUSetReuseOrdering(), PCLUSetReuseFill()
319: @*/
320: int PCILUDTSetReuseFill(PC pc,PetscTruth flag)
321: {
322:   int ierr,(*f)(PC,PetscTruth);

326:   PetscObjectQueryFunction((PetscObject)pc,"PCILUDTSetReuseFill_C",(void (**)())&f);
327:   if (f) {
328:     (*f)(pc,flag);
329:   }
330:   return(0);
331: }

333: /*@
334:    PCILUSetLevels - Sets the number of levels of fill to use.

336:    Collective on PC

338:    Input Parameters:
339: +  pc - the preconditioner context
340: -  levels - number of levels of fill

342:    Options Database Key:
343: .  -pc_ilu_levels <levels> - Sets fill level

345:    Level: intermediate

347: .keywords: PC, levels, fill, factorization, incomplete, ILU
348: @*/
349: int PCILUSetLevels(PC pc,int levels)
350: {
351:   int ierr,(*f)(PC,int);

355:   if (levels < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"negative levels");
356:   PetscObjectQueryFunction((PetscObject)pc,"PCILUSetLevels_C",(void (**)())&f);
357:   if (f) {
358:     (*f)(pc,levels);
359:   }
360:   return(0);
361: }

363: /*@
364:    PCILUSetAllowDiagonalFill - Causes all diagonal matrix entries to be 
365:    treated as level 0 fill even if there is no non-zero location.

367:    Collective on PC

369:    Input Parameters:
370: +  pc - the preconditioner context

372:    Options Database Key:
373: .  -pc_ilu_diagonal_fill

375:    Notes:
376:    Does not apply with 0 fill.

378:    Level: intermediate

380: .keywords: PC, levels, fill, factorization, incomplete, ILU
381: @*/
382: int PCILUSetAllowDiagonalFill(PC pc)
383: {
384:   int ierr,(*f)(PC);

388:   PetscObjectQueryFunction((PetscObject)pc,"PCILUSetAllowDiagonalFill_C",(void (**)())&f);
389:   if (f) {
390:     (*f)(pc);
391:   }
392:   return(0);
393: }

395: /*@
396:    PCILUSetUseInPlace - Tells the system to do an in-place incomplete factorization.
397:    Collective on PC

399:    Input Parameters:
400: .  pc - the preconditioner context

402:    Options Database Key:
403: .  -pc_ilu_in_place - Activates in-place factorization

405:    Notes:
406:    PCILUSetUseInPlace() is intended for use with matrix-free variants of
407:    Krylov methods, or when a different matrices are employed for the linear
408:    system and preconditioner, or with ASM preconditioning.  Do NOT use 
409:    this option if the linear system
410:    matrix also serves as the preconditioning matrix, since the factored
411:    matrix would then overwrite the original matrix. 

413:    Only works well with ILU(0).

415:    Level: intermediate

417: .keywords: PC, set, factorization, inplace, in-place, ILU

419: .seealso:  PCLUSetUseInPlace()
420: @*/
421: int PCILUSetUseInPlace(PC pc)
422: {
423:   int ierr,(*f)(PC);

427:   PetscObjectQueryFunction((PetscObject)pc,"PCILUSetUseInPlace_C",(void (**)())&f);
428:   if (f) {
429:     (*f)(pc);
430:   }
431:   return(0);
432: }

434: /* ------------------------------------------------------------------------------------------*/

436: static int PCSetFromOptions_ILU(PC pc)
437: {
438:   int        ierr,dtmax = 3,itmp;
439:   PetscTruth flg;
440:   PetscReal  dt[3];
441:   char       tname[256];
442:   PC_ILU     *ilu = (PC_ILU*)pc->data;
443:   PetscFList ordlist;
444:   double     tol;

447:   MatOrderingRegisterAll(PETSC_NULL);
448:   PetscOptionsHead("ILU Options");
449:     PetscOptionsInt("-pc_ilu_levels","levels of fill","PCILUSetLevels",(int)ilu->info.levels,&itmp,&flg);
450:     if (flg) ilu->info.levels = (double) itmp;
451:     PetscOptionsName("-pc_ilu_in_place","do factorization in place","PCILUSetUseInPlace",&ilu->inplace);
452:     PetscOptionsName("-pc_ilu_diagonal_fill","Allow fill into empty diagonal entry","PCILUSetAllowDiagonalFill",&flg);
453:     ilu->info.diagonal_fill = (double) flg;
454:     PetscOptionsName("-pc_iludt_reuse_fill","Reuse fill from previous ILUdt","PCILUDTSetReuseFill",&ilu->reusefill);
455:     PetscOptionsName("-pc_ilu_reuse_ordering","Reuse previous reordering","PCILUSetReuseOrdering",&ilu->reuseordering);
456:     PetscOptionsHasName(pc->prefix,"-pc_ilu_damping",&flg);
457:     if (flg) {
458:       PCILUSetDamping(pc,0.0);
459:     }
460:     PetscOptionsDouble("-pc_ilu_damping","Damping added to diagonal","PCILUSetDamping",ilu->info.damping,&ilu->info.damping,0);

462:     dt[0] = ilu->info.dt;
463:     dt[1] = ilu->info.dtcol;
464:     dt[2] = ilu->info.dtcount;
465:     PetscOptionsDoubleArray("-pc_ilu_use_drop_tolerance","<dt,dtcol,maxrowcount>","PCILUSetUseDropTolerance",dt,&dtmax,&flg);
466:     if (flg) {
467:       PCILUSetUseDropTolerance(pc,dt[0],dt[1],(int)dt[2]);
468:     }
469:     PetscOptionsDouble("-pc_ilu_fill","Expected fill in factorization","PCILUSetFill",ilu->info.fill,&ilu->info.fill,&flg);
470:     PetscOptionsDouble("-pc_ilu_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","MatReorderForNonzeroDiagonal",0.0,&tol,0);

472:     MatGetOrderingList(&ordlist);
473:     PetscOptionsList("-pc_ilu_mat_ordering_type","Reorder to reduce nonzeros in ILU","PCILUSetMatOrdering",ordlist,ilu->ordering,tname,256,&flg);
474:     if (flg) {
475:       PCILUSetMatOrdering(pc,tname);
476:     }
477:   PetscOptionsTail();
478:   return(0);
479: }

481: static int PCView_ILU(PC pc,PetscViewer viewer)
482: {
483:   PC_ILU     *ilu = (PC_ILU*)pc->data;
484:   int        ierr;
485:   PetscTruth isstring,isascii;

488:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
489:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
490:   if (isascii) {
491:     if (ilu->usedt) {
492:         PetscViewerASCIIPrintf(viewer,"  ILU: drop tolerance %gn",ilu->info.dt);
493:         PetscViewerASCIIPrintf(viewer,"  ILU: max nonzeros per row %dn",(int)ilu->info.dtcount);
494:         PetscViewerASCIIPrintf(viewer,"  ILU: column permutation tolerance %gn",ilu->info.dtcol);
495:     } else if (ilu->info.levels == 1) {
496:         PetscViewerASCIIPrintf(viewer,"  ILU: %d level of filln",(int)ilu->info.levels);
497:     } else {
498:         PetscViewerASCIIPrintf(viewer,"  ILU: %d levels of filln",(int)ilu->info.levels);
499:     }
500:     PetscViewerASCIIPrintf(viewer,"  ILU: max fill ratio allocated %gn",ilu->info.fill);
501:     if (ilu->inplace) {PetscViewerASCIIPrintf(viewer,"       in-place factorizationn");}
502:     else              {PetscViewerASCIIPrintf(viewer,"       out-of-place factorizationn");}
503:     PetscViewerASCIIPrintf(viewer,"       matrix ordering: %sn",ilu->ordering);
504:     if (ilu->reusefill)     {PetscViewerASCIIPrintf(viewer,"       Reusing fill from past factorizationn");}
505:     if (ilu->reuseordering) {PetscViewerASCIIPrintf(viewer,"       Reusing reordering from past factorizationn");}
506:     PetscViewerASCIIPrintf(viewer,"       Factored matrix followsn");
507:     PetscViewerASCIIPushTab(viewer);
508:     PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
509:     MatView(ilu->fact,viewer);
510:     PetscViewerPopFormat(viewer);
511:     PetscViewerASCIIPopTab(viewer);
512:   } else if (isstring) {
513:     PetscViewerStringSPrintf(viewer," lvls=%g,order=%s",ilu->info.levels,ilu->ordering);
514:   } else {
515:     SETERRQ1(1,"Viewer type %s not supported for PCILU",((PetscObject)viewer)->type_name);
516:   }
517:   return(0);
518: }

520: static int PCSetUp_ILU(PC pc)
521: {
522:   int        ierr;
523:   PetscTruth flg;
524:   PC_ILU     *ilu = (PC_ILU*)pc->data;

527:   if (ilu->inplace) {
528:     if (!pc->setupcalled) {

530:       /* In-place factorization only makes sense with the natural ordering,
531:          so we only need to get the ordering once, even if nonzero structure changes */
532:       MatGetOrdering(pc->pmat,ilu->ordering,&ilu->row,&ilu->col);
533:       if (ilu->row) PetscLogObjectParent(pc,ilu->row);
534:       if (ilu->col) PetscLogObjectParent(pc,ilu->col);
535:     }

537:     /* In place ILU only makes sense with fill factor of 1.0 because 
538:        cannot have levels of fill */
539:     ilu->info.fill          = 1.0;
540:     ilu->info.diagonal_fill = 0;
541:     MatILUFactor(pc->pmat,ilu->row,ilu->col,&ilu->info);
542:     ilu->fact = pc->pmat;
543:   } else if (ilu->usedt) {
544:     if (!pc->setupcalled) {
545:       MatGetOrdering(pc->pmat,ilu->ordering,&ilu->row,&ilu->col);
546:       if (ilu->row) PetscLogObjectParent(pc,ilu->row);
547:       if (ilu->col) PetscLogObjectParent(pc,ilu->col);
548:       MatILUDTFactor(pc->pmat,&ilu->info,ilu->row,ilu->col,&ilu->fact);
549:       PetscLogObjectParent(pc,ilu->fact);
550:     } else if (pc->flag != SAME_NONZERO_PATTERN) {
551:       MatDestroy(ilu->fact);
552:       if (!ilu->reuseordering) {
553:         if (ilu->row) {ISDestroy(ilu->row);}
554:         if (ilu->col) {ISDestroy(ilu->col);}
555:         MatGetOrdering(pc->pmat,ilu->ordering,&ilu->row,&ilu->col);
556:         if (ilu->row) PetscLogObjectParent(pc,ilu->row);
557:         if (ilu->col) PetscLogObjectParent(pc,ilu->col);
558:       }
559:       MatILUDTFactor(pc->pmat,&ilu->info,ilu->row,ilu->col,&ilu->fact);
560:       PetscLogObjectParent(pc,ilu->fact);
561:     } else if (!ilu->reusefill) {
562:       MatDestroy(ilu->fact);
563:       MatILUDTFactor(pc->pmat,&ilu->info,ilu->row,ilu->col,&ilu->fact);
564:       PetscLogObjectParent(pc,ilu->fact);
565:     } else {
566:       MatLUFactorNumeric(pc->pmat,&ilu->fact);
567:     }
568:   } else {
569:     if (!pc->setupcalled) {
570:       /* first time in so compute reordering and symbolic factorization */
571:       MatGetOrdering(pc->pmat,ilu->ordering,&ilu->row,&ilu->col);
572:       if (ilu->row) PetscLogObjectParent(pc,ilu->row);
573:       if (ilu->col) PetscLogObjectParent(pc,ilu->col);
574:       /*  Remove zeros along diagonal?     */
575:       PetscOptionsHasName(pc->prefix,"-pc_ilu_nonzeros_along_diagonal",&flg);
576:       if (flg) {
577:         PetscReal ntol = 1.e-10;
578:         PetscOptionsGetDouble(pc->prefix,"-pc_ilu_nonzeros_along_diagonal",&ntol,PETSC_NULL);
579:         MatReorderForNonzeroDiagonal(pc->pmat,ntol,ilu->row,ilu->col);
580:       }

582:       MatILUFactorSymbolic(pc->pmat,ilu->row,ilu->col,&ilu->info,&ilu->fact);
583:       PetscLogObjectParent(pc,ilu->fact);
584:     } else if (pc->flag != SAME_NONZERO_PATTERN) {
585:       if (!ilu->reuseordering) {
586:         /* compute a new ordering for the ILU */
587:         ISDestroy(ilu->row);
588:         ISDestroy(ilu->col);
589:         MatGetOrdering(pc->pmat,ilu->ordering,&ilu->row,&ilu->col);
590:         if (ilu->row) PetscLogObjectParent(pc,ilu->row);
591:         if (ilu->col) PetscLogObjectParent(pc,ilu->col);
592:         /*  Remove zeros along diagonal?     */
593:         PetscOptionsHasName(pc->prefix,"-pc_ilu_nonzeros_along_diagonal",&flg);
594:         if (flg) {
595:           PetscReal ntol = 1.e-10;
596:           PetscOptionsGetDouble(pc->prefix,"-pc_ilu_nonzeros_along_diagonal",&ntol,PETSC_NULL);
597:           MatReorderForNonzeroDiagonal(pc->pmat,ntol,ilu->row,ilu->col);
598:         }
599:       }
600:       MatDestroy(ilu->fact);
601:       MatILUFactorSymbolic(pc->pmat,ilu->row,ilu->col,&ilu->info,&ilu->fact);
602:       PetscLogObjectParent(pc,ilu->fact);
603:     }
604:     MatLUFactorNumeric(pc->pmat,&ilu->fact);
605:   }
606:   return(0);
607: }

609: static int PCDestroy_ILU(PC pc)
610: {
611:   PC_ILU *ilu = (PC_ILU*)pc->data;
612:   int    ierr;

615:   if (!ilu->inplace && ilu->fact) {MatDestroy(ilu->fact);}
616:   if (ilu->row && ilu->col && ilu->row != ilu->col) {ISDestroy(ilu->row);}
617:   if (ilu->col) {ISDestroy(ilu->col);}
618:   PetscStrfree(ilu->ordering);
619:   PetscFree(ilu);
620:   return(0);
621: }

623: static int PCApply_ILU(PC pc,Vec x,Vec y)
624: {
625:   PC_ILU *ilu = (PC_ILU*)pc->data;
626:   int    ierr;

629:   MatSolve(ilu->fact,x,y);
630:   return(0);
631: }

633: static int PCApplyTranspose_ILU(PC pc,Vec x,Vec y)
634: {
635:   PC_ILU *ilu = (PC_ILU*)pc->data;
636:   int    ierr;

639:   MatSolveTranspose(ilu->fact,x,y);
640:   return(0);
641: }

643: static int PCGetFactoredMatrix_ILU(PC pc,Mat *mat)
644: {
645:   PC_ILU *ilu = (PC_ILU*)pc->data;

648:   if (!ilu->fact) SETERRQ(1,"Matrix not yet factored; call after SLESSetUp() or PCSetUp()");
649:   *mat = ilu->fact;
650:   return(0);
651: }

653: EXTERN_C_BEGIN
654: int PCCreate_ILU(PC pc)
655: {
656:   int    ierr;
657:   PC_ILU *ilu;

660:   PetscNew(PC_ILU,&ilu);
661:   PetscLogObjectMemory(pc,sizeof(PC_ILU));

663:   ilu->fact               = 0;
664:   ilu->info.levels        = 0;
665:   ilu->info.fill          = 1.0;
666:   ilu->col                = 0;
667:   ilu->row                = 0;
668:   ilu->inplace            = PETSC_FALSE;
669:   PetscStrallocpy(MATORDERING_NATURAL,&ilu->ordering);
670:   ilu->reuseordering      = PETSC_FALSE;
671:   ilu->usedt              = PETSC_FALSE;
672:   ilu->info.dt            = PETSC_DEFAULT;
673:   ilu->info.dtcount       = PETSC_DEFAULT;
674:   ilu->info.dtcol         = PETSC_DEFAULT;
675:   ilu->info.damp          = 0.0;
676:   ilu->info.damping       = 0.0;
677:   ilu->reusefill          = PETSC_FALSE;
678:   ilu->info.diagonal_fill = 0;
679:   pc->data                = (void*)ilu;

681:   pc->ops->destroy           = PCDestroy_ILU;
682:   pc->ops->apply             = PCApply_ILU;
683:   pc->ops->applytranspose    = PCApplyTranspose_ILU;
684:   pc->ops->setup             = PCSetUp_ILU;
685:   pc->ops->setfromoptions    = PCSetFromOptions_ILU;
686:   pc->ops->getfactoredmatrix = PCGetFactoredMatrix_ILU;
687:   pc->ops->view              = PCView_ILU;
688:   pc->ops->applyrichardson   = 0;

690:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCILUSetUseDropTolerance_C","PCILUSetUseDropTolerance_ILU",
691:                     PCILUSetUseDropTolerance_ILU);
692:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCILUSetFill_C","PCILUSetFill_ILU",
693:                     PCILUSetFill_ILU);
694:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCILUSetDamping_C","PCILUSetDamping_ILU",
695:                     PCILUSetDamping_ILU);
696:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCILUSetMatOrdering_C","PCILUSetMatOrdering_ILU",
697:                     PCILUSetMatOrdering_ILU);
698:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCILUSetReuseOrdering_C","PCILUSetReuseOrdering_ILU",
699:                     PCILUSetReuseOrdering_ILU);
700:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCILUDTSetReuseFill_C","PCILUDTSetReuseFill_ILUDT",
701:                     PCILUDTSetReuseFill_ILUDT);
702:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCILUSetLevels_C","PCILUSetLevels_ILU",
703:                     PCILUSetLevels_ILU);
704:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCILUSetUseInPlace_C","PCILUSetUseInPlace_ILU",
705:                     PCILUSetUseInPlace_ILU);
706:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCILUSetAllowDiagonalFill_C","PCILUSetAllowDiagonalFill_ILU",
707:                     PCILUSetAllowDiagonalFill_ILU);

709:   return(0);
710: }
711: EXTERN_C_END