Actual source code: lu.c

  1: /*$Id: lu.c,v 1.149 2001/08/07 21:30:26 bsmith Exp $*/
  2: /*
  3:    Defines a direct factorization preconditioner for any Mat implementation
  4:    Note: this need not be consided a preconditioner since it supplies
  5:          a direct solver.
  6: */
 7:  #include src/sles/pc/pcimpl.h

  9: typedef struct {
 10:   Mat             fact;             /* factored matrix */
 11:   PetscReal       actualfill;       /* actual fill in factor */
 12:   int             inplace;          /* flag indicating in-place factorization */
 13:   IS              row,col;          /* index sets used for reordering */
 14:   MatOrderingType ordering;         /* matrix ordering */
 15:   PetscTruth      reuseordering;    /* reuses previous reordering computed */
 16:   PetscTruth      reusefill;        /* reuse fill from previous LU */
 17:   MatLUInfo       info;
 18: } PC_LU;


 21: EXTERN_C_BEGIN
 22: #undef __FUNCT__  
 24: int PCLUSetZeroPivot_LU(PC pc,PetscReal z)
 25: {
 26:   PC_LU *lu;

 29:   lu                 = (PC_LU*)pc->data;
 30:   lu->info.zeropivot = z;
 31:   return(0);
 32: }
 33: EXTERN_C_END

 35: EXTERN_C_BEGIN
 36: #undef __FUNCT__  
 38: int PCLUSetReuseOrdering_LU(PC pc,PetscTruth flag)
 39: {
 40:   PC_LU *lu;

 43:   lu                = (PC_LU*)pc->data;
 44:   lu->reuseordering = flag;
 45:   return(0);
 46: }
 47: EXTERN_C_END

 49: EXTERN_C_BEGIN
 50: #undef __FUNCT__  
 52: int PCLUSetReuseFill_LU(PC pc,PetscTruth flag)
 53: {
 54:   PC_LU *lu;

 57:   lu = (PC_LU*)pc->data;
 58:   lu->reusefill = flag;
 59:   return(0);
 60: }
 61: EXTERN_C_END

 63: #undef __FUNCT__  
 65: static int PCSetFromOptions_LU(PC pc)
 66: {
 67:   PC_LU      *lu = (PC_LU*)pc->data;
 68:   int        ierr;
 69:   PetscTruth flg,set;
 70:   char       tname[256];
 71:   PetscFList ordlist;
 72:   PetscReal  tol;

 75:   MatOrderingRegisterAll(PETSC_NULL);
 76:   PetscOptionsHead("LU options");
 77:     PetscOptionsName("-pc_lu_in_place","Form LU in the same memory as the matrix","PCLUSetUseInPlace",&flg);
 78:     if (flg) {
 79:       PCLUSetUseInPlace(pc);
 80:     }
 81:     PetscOptionsReal("-pc_lu_fill","Expected non-zeros in LU/non-zeros in matrix","PCLUSetFill",lu->info.fill,&lu->info.fill,0);

 83:     PetscOptionsName("-pc_lu_damping","Damping added to diagonal","PCLUSetDamping",&flg);
 84:     if (flg) {
 85:         PCLUSetDamping(pc,(PetscReal) PETSC_DECIDE);
 86:     }
 87:     PetscOptionsReal("-pc_lu_damping","Damping added to diagonal","PCLUSetDamping",lu->info.damping,&lu->info.damping,0);
 88:     PetscOptionsReal("-pc_lu_zeropivot","Pivot is considered zero if less than","PCLUSetSetZeroPivot",lu->info.zeropivot,&lu->info.zeropivot,0);

 90:     PetscOptionsName("-pc_lu_reuse_fill","Use fill from previous factorization","PCLUSetReuseFill",&flg);
 91:     if (flg) {
 92:       PCLUSetReuseFill(pc,PETSC_TRUE);
 93:     }
 94:     PetscOptionsName("-pc_lu_reuse_ordering","Reuse ordering from previous factorization","PCLUSetReuseOrdering",&flg);
 95:     if (flg) {
 96:       PCLUSetReuseOrdering(pc,PETSC_TRUE);
 97:     }

 99:     MatGetOrderingList(&ordlist);
100:     PetscOptionsList("-pc_lu_mat_ordering_type","Reordering to reduce nonzeros in LU","PCLUSetMatOrdering",ordlist,lu->ordering,tname,256,&flg);
101:     if (flg) {
102:       PCLUSetMatOrdering(pc,tname);
103:     }
104:     PetscOptionsReal("-pc_lu_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","MatReorderForNonzeroDiagonal",0.0,&tol,0);

106:     PetscOptionsReal("-pc_lu_pivoting","Pivoting tolerance (used only for some factorization)","PCLUSetPivoting",lu->info.dtcol,&lu->info.dtcol,&flg);

108:     flg = lu->info.pivotinblocks ? PETSC_TRUE : PETSC_FALSE;
109:     PetscOptionsLogical("-pc_lu_pivot_in_blocks","Pivot inside matrix blocks for BAIJ and SBAIJ","PCLUSetPivotInBlocks",flg,&flg,&set);
110:     if (set) {
111:       PCLUSetPivotInBlocks(pc,flg);
112:     }
113:   PetscOptionsTail();
114:   return(0);
115: }

117: #undef __FUNCT__  
119: static int PCView_LU(PC pc,PetscViewer viewer)
120: {
121:   PC_LU      *lu = (PC_LU*)pc->data;
122:   int        ierr;
123:   PetscTruth isascii,isstring;

126:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
127:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
128:   if (isascii) {
129:     MatInfo info;

131:     if (lu->inplace) {PetscViewerASCIIPrintf(viewer,"  LU: in-place factorizationn");}
132:     else             {PetscViewerASCIIPrintf(viewer,"  LU: out-of-place factorizationn");}
133:     PetscViewerASCIIPrintf(viewer,"    matrix ordering: %sn",lu->ordering);
134:     PetscViewerASCIIPrintf(viewer,"  LU: tolerance for zero pivot %gn",lu->info.zeropivot);
135:     if (lu->fact) {
136:       MatGetInfo(lu->fact,MAT_LOCAL,&info);
137:       PetscViewerASCIIPrintf(viewer,"    LU nonzeros %gn",info.nz_used);
138:       PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_FACTOR_INFO);
139:       MatView(lu->fact,viewer);
140:       PetscViewerPopFormat(viewer);
141:     }
142:     if (lu->reusefill)    {PetscViewerASCIIPrintf(viewer,"       Reusing fill from past factorizationn");}
143:     if (lu->reuseordering) {PetscViewerASCIIPrintf(viewer,"       Reusing reordering from past factorizationn");}
144:   } else if (isstring) {
145:     PetscViewerStringSPrintf(viewer," order=%s",lu->ordering);
146:   } else {
147:     SETERRQ1(1,"Viewer type %s not supported for PCLU",((PetscObject)viewer)->type_name);
148:   }
149:   return(0);
150: }

152: #undef __FUNCT__  
154: static int PCGetFactoredMatrix_LU(PC pc,Mat *mat)
155: {
156:   PC_LU *dir = (PC_LU*)pc->data;

159:   if (!dir->fact) SETERRQ(1,"Matrix not yet factored; call after SLESSetUp() or PCSetUp()");
160:   *mat = dir->fact;
161:   return(0);
162: }

164: #undef __FUNCT__  
166: static int PCSetUp_LU(PC pc)
167: {
168:   int        ierr;
169:   PetscTruth flg;
170:   PC_LU      *dir = (PC_LU*)pc->data;

173:   if (dir->reusefill && pc->setupcalled) dir->info.fill = dir->actualfill;

175:   if (dir->inplace) {
176:     if (dir->row && dir->col && dir->row != dir->col) {ISDestroy(dir->row);}
177:     if (dir->col) {ISDestroy(dir->col);}
178:     MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);
179:     if (dir->row) {PetscLogObjectParent(pc,dir->row); PetscLogObjectParent(pc,dir->col);}
180:     MatLUFactor(pc->pmat,dir->row,dir->col,&dir->info);
181:     dir->fact = pc->pmat;
182:   } else {
183:     MatInfo info;
184:     if (!pc->setupcalled) {
185:       MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);
186:       PetscOptionsHasName(pc->prefix,"-pc_lu_nonzeros_along_diagonal",&flg);
187:       if (flg) {
188:         PetscReal tol = 1.e-10;
189:         PetscOptionsGetReal(pc->prefix,"-pc_lu_nonzeros_along_diagonal",&tol,PETSC_NULL);
190:         MatReorderForNonzeroDiagonal(pc->pmat,tol,dir->row,dir->col);
191:       }
192:       if (dir->row) {PetscLogObjectParent(pc,dir->row); PetscLogObjectParent(pc,dir->col);}
193:       MatLUFactorSymbolic(pc->pmat,dir->row,dir->col,&dir->info,&dir->fact);
194:       MatGetInfo(dir->fact,MAT_LOCAL,&info);
195:       dir->actualfill = info.fill_ratio_needed;
196:       PetscLogObjectParent(pc,dir->fact);
197:     } else if (pc->flag != SAME_NONZERO_PATTERN) {
198:       if (!dir->reuseordering) {
199:         if (dir->row && dir->col && dir->row != dir->col) {ISDestroy(dir->row);}
200:         if (dir->col) {ISDestroy(dir->col);}
201:         MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);
202:         PetscOptionsHasName(pc->prefix,"-pc_lu_nonzeros_along_diagonal",&flg);
203:         if (flg) {
204:           PetscReal tol = 1.e-10;
205:           PetscOptionsGetReal(pc->prefix,"-pc_lu_nonzeros_along_diagonal",&tol,PETSC_NULL);
206:           MatReorderForNonzeroDiagonal(pc->pmat,tol,dir->row,dir->col);
207:         }
208:         if (dir->row) {PetscLogObjectParent(pc,dir->row); PetscLogObjectParent(pc,dir->col);}
209:       }
210:       MatDestroy(dir->fact);
211:       MatLUFactorSymbolic(pc->pmat,dir->row,dir->col,&dir->info,&dir->fact);
212:       MatGetInfo(dir->fact,MAT_LOCAL,&info);
213:       dir->actualfill = info.fill_ratio_needed;
214:       PetscLogObjectParent(pc,dir->fact);
215:     }
216:     MatLUFactorNumeric(pc->pmat,&dir->fact);
217:   }
218:   return(0);
219: }

221: #undef __FUNCT__  
223: static int PCDestroy_LU(PC pc)
224: {
225:   PC_LU *dir = (PC_LU*)pc->data;
226:   int   ierr;

229:   if (!dir->inplace && dir->fact) {MatDestroy(dir->fact);}
230:   if (dir->row && dir->col && dir->row != dir->col) {ISDestroy(dir->row);}
231:   if (dir->col) {ISDestroy(dir->col);}
232:   PetscStrfree(dir->ordering);
233:   PetscFree(dir);
234:   return(0);
235: }

237: #undef __FUNCT__  
239: static int PCApply_LU(PC pc,Vec x,Vec y)
240: {
241:   PC_LU *dir = (PC_LU*)pc->data;
242:   int   ierr;

245:   if (dir->inplace) {MatSolve(pc->pmat,x,y);}
246:   else              {MatSolve(dir->fact,x,y);}
247:   return(0);
248: }

250: #undef __FUNCT__  
252: static int PCApplyTranspose_LU(PC pc,Vec x,Vec y)
253: {
254:   PC_LU *dir = (PC_LU*)pc->data;
255:   int   ierr;

258:   if (dir->inplace) {MatSolveTranspose(pc->pmat,x,y);}
259:   else              {MatSolveTranspose(dir->fact,x,y);}
260:   return(0);
261: }

263: /* -----------------------------------------------------------------------------------*/

265: EXTERN_C_BEGIN
266: #undef __FUNCT__  
268: int PCLUSetFill_LU(PC pc,PetscReal fill)
269: {
270:   PC_LU *dir;

273:   dir = (PC_LU*)pc->data;
274:   dir->info.fill = fill;
275:   return(0);
276: }
277: EXTERN_C_END

279: EXTERN_C_BEGIN
280: #undef __FUNCT__  
282: int PCLUSetDamping_LU(PC pc,PetscReal damping)
283: {
284:   PC_LU *dir;

287:   dir = (PC_LU*)pc->data;
288:   if (damping == (PetscReal) PETSC_DECIDE) {
289:     dir->info.damping = 1.e-12;
290:   } else {
291:     dir->info.damping = damping;
292:   }
293:   return(0);
294: }
295: EXTERN_C_END

297: EXTERN_C_BEGIN
298: #undef __FUNCT__  
300: int PCLUSetUseInPlace_LU(PC pc)
301: {
302:   PC_LU *dir;

305:   dir = (PC_LU*)pc->data;
306:   dir->inplace = 1;
307:   return(0);
308: }
309: EXTERN_C_END

311: EXTERN_C_BEGIN
312: #undef __FUNCT__  
314: int PCLUSetMatOrdering_LU(PC pc,MatOrderingType ordering)
315: {
316:   PC_LU *dir = (PC_LU*)pc->data;
317:   int   ierr;

320:   PetscStrfree(dir->ordering);
321:   PetscStrallocpy(ordering,&dir->ordering);
322:   return(0);
323: }
324: EXTERN_C_END

326: EXTERN_C_BEGIN
327: #undef __FUNCT__  
329: int PCLUSetPivoting_LU(PC pc,PetscReal dtcol)
330: {
331:   PC_LU *dir = (PC_LU*)pc->data;

334:   if (dtcol < 0.0 || dtcol > 1.0) SETERRQ1(1,"Column pivot tolerance is %g must be between 0 and 1",dtcol);
335:   dir->info.dtcol = dtcol;
336:   return(0);
337: }
338: EXTERN_C_END

340: EXTERN_C_BEGIN
341: #undef __FUNCT__  
343: int PCLUSetPivotInBlocks_LU(PC pc,PetscTruth pivot)
344: {
345:   PC_LU *dir = (PC_LU*)pc->data;

348:   dir->info.pivotinblocks = pivot ? 1.0 : 0.0;
349:   return(0);
350: }
351: EXTERN_C_END

353: /* -----------------------------------------------------------------------------------*/

355: #undef __FUNCT__  
357: /*@
358:    PCLUSetZeroPivot - Sets the size at which smaller pivots are declared to be zero

360:    Collective on PC
361:    
362:    Input Parameters:
363: +  pc - the preconditioner context
364: -  zero - all pivots smaller than this will be considered zero

366:    Options Database Key:
367: .  -pc_ilu_zeropivot <zero> - Sets the zero pivot size

369:    Level: intermediate

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

373: .seealso: PCLUSetFill(), PCLUSetDamp(), PCILUSetZeroPivot()
374: @*/
375: int PCLUSetZeroPivot(PC pc,PetscReal zero)
376: {
377:   int ierr,(*f)(PC,PetscReal);

381:   PetscObjectQueryFunction((PetscObject)pc,"PCLUSetZeroPivot_C",(void (**)(void))&f);
382:   if (f) {
383:     (*f)(pc,zero);
384:   }
385:   return(0);
386: }

388: #undef __FUNCT__  
390: /*@
391:    PCLUSetReuseOrdering - When similar matrices are factored, this
392:    causes the ordering computed in the first factor to be used for all
393:    following factors; applies to both fill and drop tolerance LUs.

395:    Collective on PC

397:    Input Parameters:
398: +  pc - the preconditioner context
399: -  flag - PETSC_TRUE to reuse else PETSC_FALSE

401:    Options Database Key:
402: .  -pc_lu_reuse_ordering - Activate PCLUSetReuseOrdering()

404:    Level: intermediate

406: .keywords: PC, levels, reordering, factorization, incomplete, LU

408: .seealso: PCLUSetReuseFill(), PCILUSetReuseOrdering(), PCILUDTSetReuseFill()
409: @*/
410: int PCLUSetReuseOrdering(PC pc,PetscTruth flag)
411: {
412:   int ierr,(*f)(PC,PetscTruth);

416:   PetscObjectQueryFunction((PetscObject)pc,"PCLUSetReuseOrdering_C",(void (**)(void))&f);
417:   if (f) {
418:     (*f)(pc,flag);
419:   }
420:   return(0);
421: }

423: #undef __FUNCT__  
425: /*@
426:    PCLUSetReuseFill - When matrices with same nonzero structure are LU factored,
427:    this causes later ones to use the fill computed in the initial factorization.

429:    Collective on PC

431:    Input Parameters:
432: +  pc - the preconditioner context
433: -  flag - PETSC_TRUE to reuse else PETSC_FALSE

435:    Options Database Key:
436: .  -pc_lu_reuse_fill - Activates PCLUSetReuseFill()

438:    Level: intermediate

440: .keywords: PC, levels, reordering, factorization, incomplete, LU

442: .seealso: PCILUSetReuseOrdering(), PCLUSetReuseOrdering(), PCILUDTSetReuseFill()
443: @*/
444: int PCLUSetReuseFill(PC pc,PetscTruth flag)
445: {
446:   int ierr,(*f)(PC,PetscTruth);

450:   PetscObjectQueryFunction((PetscObject)pc,"PCLUSetReuseFill_C",(void (**)(void))&f);
451:   if (f) {
452:     (*f)(pc,flag);
453:   }
454:   return(0);
455: }

457: #undef __FUNCT__  
459: /*@
460:    PCLUSetFill - Indicate the amount of fill you expect in the factored matrix,
461:    fill = number nonzeros in factor/number nonzeros in original matrix.

463:    Collective on PC
464:    
465:    Input Parameters:
466: +  pc - the preconditioner context
467: -  fill - amount of expected fill

469:    Options Database Key:
470: .  -pc_lu_fill <fill> - Sets fill amount

472:    Level: intermediate

474:    Note:
475:    For sparse matrix factorizations it is difficult to predict how much 
476:    fill to expect. By running with the option -log_info PETSc will print the 
477:    actual amount of fill used; allowing you to set the value accurately for
478:    future runs. Default PETSc uses a value of 5.0

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

482: .seealso: PCILUSetFill()
483: @*/
484: int PCLUSetFill(PC pc,PetscReal fill)
485: {
486:   int ierr,(*f)(PC,PetscReal);

490:   if (fill < 1.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Fill factor cannot be less then 1.0");
491:   PetscObjectQueryFunction((PetscObject)pc,"PCLUSetFill_C",(void (**)(void))&f);
492:   if (f) {
493:     (*f)(pc,fill);
494:   }
495:   return(0);
496: }

498: #undef __FUNCT__  
500: /*@
501:    PCLUSetDamping - adds this quantity to the diagonal of the matrix during the 
502:      LU numerical factorization

504:    Collective on PC
505:    
506:    Input Parameters:
507: +  pc - the preconditioner context
508: -  damping - amount of damping  (use PETSC_DECIDE for default of 1.e-12)

510:    Options Database Key:
511: .  -pc_lu_damping <damping> - Sets damping amount or PETSC_DECIDE for the default

513:    Note: If 0.0 is given, then no damping is used. If a diagonal element is classified as a zero
514:          pivot, then the damping is doubled until this is alleviated.

516:    Level: intermediate

518: .keywords: PC, set, factorization, direct, fill
519: .seealso: PCILUSetFill(), PCILUSetDamp()
520: @*/
521: int PCLUSetDamping(PC pc,PetscReal damping)
522: {
523:   int ierr,(*f)(PC,PetscReal);

527:   PetscObjectQueryFunction((PetscObject)pc,"PCLUSetDamping_C",(void (**)(void))&f);
528:   if (f) {
529:     (*f)(pc,damping);
530:   }
531:   return(0);
532: }

534: #undef __FUNCT__  
536: /*@
537:    PCLUSetUseInPlace - Tells the system to do an in-place factorization.
538:    For dense matrices, this enables the solution of much larger problems. 
539:    For sparse matrices the factorization cannot be done truly in-place 
540:    so this does not save memory during the factorization, but after the matrix
541:    is factored, the original unfactored matrix is freed, thus recovering that
542:    space.

544:    Collective on PC

546:    Input Parameters:
547: .  pc - the preconditioner context

549:    Options Database Key:
550: .  -pc_lu_in_place - Activates in-place factorization

552:    Notes:
553:    PCLUSetUseInplace() can only be used with the KSP method KSPPREONLY or when 
554:    a different matrix is provided for the multiply and the preconditioner in 
555:    a call to SLESSetOperators().
556:    This is because the Krylov space methods require an application of the 
557:    matrix multiplication, which is not possible here because the matrix has 
558:    been factored in-place, replacing the original matrix.

560:    Level: intermediate

562: .keywords: PC, set, factorization, direct, inplace, in-place, LU

564: .seealso: PCILUSetUseInPlace()
565: @*/
566: int PCLUSetUseInPlace(PC pc)
567: {
568:   int ierr,(*f)(PC);

572:   PetscObjectQueryFunction((PetscObject)pc,"PCLUSetUseInPlace_C",(void (**)(void))&f);
573:   if (f) {
574:     (*f)(pc);
575:   }
576:   return(0);
577: }

579: #undef __FUNCT__  
581: /*@C
582:     PCLUSetMatOrdering - Sets the ordering routine (to reduce fill) to 
583:     be used in the LU factorization.

585:     Collective on PC

587:     Input Parameters:
588: +   pc - the preconditioner context
589: -   ordering - the matrix ordering name, for example, MATORDERING_ND or MATORDERING_RCM

591:     Options Database Key:
592: .   -pc_lu_mat_ordering_type <nd,rcm,...> - Sets ordering routine

594:     Level: intermediate

596:     Notes: nested dissection is used by default

598: .seealso: PCILUSetMatOrdering()
599: @*/
600: int PCLUSetMatOrdering(PC pc,MatOrderingType ordering)
601: {
602:   int ierr,(*f)(PC,MatOrderingType);

605:   PetscObjectQueryFunction((PetscObject)pc,"PCLUSetMatOrdering_C",(void (**)(void))&f);
606:   if (f) {
607:     (*f)(pc,ordering);
608:   }
609:   return(0);
610: }

612: #undef __FUNCT__  
614: /*@
615:     PCLUSetPivoting - Determines when pivoting is done during LU. 
616:       For PETSc dense matrices column pivoting is always done, for PETSc sparse matrices
617:       it is never done. For the Matlab and SuperLU factorization this is used.

619:     Collective on PC

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

625:     Options Database Key:
626: .   -pc_lu_pivoting - dttol

628:     Level: intermediate

630: .seealso: PCILUSetMatOrdering(), PCLUSetPivotInBlocks()
631: @*/
632: int PCLUSetPivoting(PC pc,PetscReal dtcol)
633: {
634:   int ierr,(*f)(PC,PetscReal);

637:   PetscObjectQueryFunction((PetscObject)pc,"PCLUSetPivoting_C",(void (**)(void))&f);
638:   if (f) {
639:     (*f)(pc,dtcol);
640:   }
641:   return(0);
642: }

644: #undef __FUNCT__  
646: /*@
647:     PCLUSetPivotInBlocks - Determines if pivoting is done while factoring each block
648:       with BAIJ or SBAIJ matrices

650:     Collective on PC

652:     Input Parameters:
653: +   pc - the preconditioner context
654: -   pivot - PETSC_TRUE or PETSC_FALSE

656:     Options Database Key:
657: .   -pc_lu_pivot_in_blocks <true,false>

659:     Level: intermediate

661: .seealso: PCILUSetMatOrdering(), PCLUSetPivoting()
662: @*/
663: int PCLUSetPivotInBlocks(PC pc,PetscTruth pivot)
664: {
665:   int ierr,(*f)(PC,PetscTruth);

668:   PetscObjectQueryFunction((PetscObject)pc,"PCLUSetPivotInBlocks_C",(void (**)(void))&f);
669:   if (f) {
670:     (*f)(pc,pivot);
671:   }
672:   return(0);
673: }

675: /* ------------------------------------------------------------------------ */

677: /*MC
678:    PCLU - Uses a direct solver, based on LU factorization, as a preconditioner

680:    Options Database Keys:
681: +  -pc_lu_reuse_ordering - Activate PCLUSetReuseOrdering()
682: .  -pc_lu_reuse_fill - Activates PCLUSetReuseFill()
683: .  -pc_lu_fill <fill> - Sets fill amount
684: .  -pc_lu_damping <damping> - Sets damping amount
685: .  -pc_lu_in_place - Activates in-place factorization
686: .  -pc_lu_mat_ordering_type <nd,rcm,...> - Sets ordering routine
687: .  -pc_lu_pivoting - dttol
688: -  -pc_lu_pivot_in_blocks <true,false>

690:    Notes: Not all options work for all matrix formats
691:           Run with -help to see additional options for particular matrix formats or factorization
692:           algorithms

694:    Level: beginner

696:    Concepts: LU factorization, direct solver

698:    Notes: Usually this will compute an "exact" solution in one iteration and does 
699:           not need a Krylov method (i.e. you can use -ksp_type preonly, or 
700:           KSPSetType(ksp,KSPPREONLY) for the Krylov method

702: .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PC,
703:            PCILU, PCCHOLESKY, PCICC, PCLUSetReuseOrdering(), PCLUSetReuseFill(), PCGetFactoredMatrix(),
704:            PCLUSetFill(), PCLUSetDamping(), PCLUSetUseInPlace(), PCLUSetMatOrdering(), PCLUSetPivoting(),
705:            PCLUSetPivotingInBlocks()
706: M*/

708: EXTERN_C_BEGIN
709: #undef __FUNCT__  
711: int PCCreate_LU(PC pc)
712: {
713:   int   ierr,size;
714:   PC_LU *dir;

717:   PetscNew(PC_LU,&dir);
718:   PetscLogObjectMemory(pc,sizeof(PC_LU));

720:   dir->fact               = 0;
721:   dir->inplace            = 0;
722:   dir->info.fill          = 5.0;
723:   dir->info.dtcol         = 1.e-6; /* default to pivoting; this is only thing PETSc LU supports */
724:   dir->info.damping       = 0.0;
725:   dir->info.zeropivot     = 1.e-12;
726:   dir->info.pivotinblocks = 1.0;
727:   dir->col                = 0;
728:   dir->row                = 0;
729:   MPI_Comm_size(pc->comm,&size);
730:   if (size == 1) {
731:     PetscStrallocpy(MATORDERING_ND,&dir->ordering);
732:   } else {
733:     PetscStrallocpy(MATORDERING_NATURAL,&dir->ordering);
734:   }
735:   dir->reusefill        = PETSC_FALSE;
736:   dir->reuseordering    = PETSC_FALSE;
737:   pc->data              = (void*)dir;

739:   pc->ops->destroy           = PCDestroy_LU;
740:   pc->ops->apply             = PCApply_LU;
741:   pc->ops->applytranspose    = PCApplyTranspose_LU;
742:   pc->ops->setup             = PCSetUp_LU;
743:   pc->ops->setfromoptions    = PCSetFromOptions_LU;
744:   pc->ops->view              = PCView_LU;
745:   pc->ops->applyrichardson   = 0;
746:   pc->ops->getfactoredmatrix = PCGetFactoredMatrix_LU;

748:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCLUSetFill_C","PCLUSetFill_LU",
749:                     PCLUSetFill_LU);
750:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCLUSetDamping_C","PCLUSetDamping_LU",
751:                     PCLUSetDamping_LU);
752:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCLUSetUseInPlace_C","PCLUSetUseInPlace_LU",
753:                     PCLUSetUseInPlace_LU);
754:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCLUSetMatOrdering_C","PCLUSetMatOrdering_LU",
755:                     PCLUSetMatOrdering_LU);
756:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCLUSetReuseOrdering_C","PCLUSetReuseOrdering_LU",
757:                     PCLUSetReuseOrdering_LU);
758:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCLUSetReuseFill_C","PCLUSetReuseFill_LU",
759:                     PCLUSetReuseFill_LU);
760:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCLUSetPivoting_C","PCLUSetPivoting_LU",
761:                     PCLUSetPivoting_LU);
762:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCLUSetPivotInBlocks_C","PCLUSetPivotInBlocks_LU",
763:                     PCLUSetPivotInBlocks_LU);
764:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCLUSetZeroPivot_C","PCLUSetZeroPivot_LU",
765:                     PCLUSetZeroPivot_LU);
766:   return(0);
767: }
768: EXTERN_C_END