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