Actual source code: icc.c
1: /*$Id: icc.c,v 1.82 2001/08/06 21:16:31 bsmith Exp $*/
2: /*
3: Defines a Cholesky factorization preconditioner for any Mat implementation.
4: Presently only provided for MPIRowbs format (i.e. BlockSolve).
5: */
7: #include src/sles/pc/impls/icc/icc.h
9: EXTERN_C_BEGIN
10: int PCICCSetMatOrdering_ICC(PC pc,MatOrderingType ordering)
11: {
12: PC_ICC *dir = (PC_ICC*)pc->data;
13: int ierr;
14:
16: PetscStrfree(dir->ordering);
17: PetscStrallocpy(ordering,&dir->ordering);
18: return(0);
19: }
20: EXTERN_C_END
22: EXTERN_C_BEGIN
23: int PCICCSetFill_ICC(PC pc,PetscReal fill)
24: {
25: PC_ICC *dir;
28: dir = (PC_ICC*)pc->data;
29: dir->fill = fill;
30: return(0);
31: }
32: EXTERN_C_END
34: EXTERN_C_BEGIN
35: int PCICCSetLevels_ICC(PC pc,int levels)
36: {
37: PC_ICC *icc;
40: icc = (PC_ICC*)pc->data;
41: icc->levels = levels;
42: return(0);
43: }
44: EXTERN_C_END
46: /*@
47: PCICCSetMatOrdering - Sets the ordering routine (to reduce fill) to
48: be used it the ICC factorization.
50: Collective on PC
52: Input Parameters:
53: + pc - the preconditioner context
54: - ordering - the matrix ordering name, for example, MATORDERING_ND or MATORDERING_RCM
56: Options Database Key:
57: . -pc_icc_mat_ordering_type <nd,rcm,...> - Sets ordering routine
59: Level: intermediate
61: .seealso: PCLUSetMatOrdering()
63: .keywords: PC, ICC, set, matrix, reordering
65: @*/
66: int PCICCSetMatOrdering(PC pc,MatOrderingType ordering)
67: {
68: int ierr,(*f)(PC,MatOrderingType);
72: PetscObjectQueryFunction((PetscObject)pc,"PCICCSetMatOrdering_C",(void (**)(void))&f);
73: if (f) {
74: (*f)(pc,ordering);
75: }
76: return(0);
77: }
79: /*@
80: PCICCSetLevels - Sets the number of levels of fill to use.
82: Collective on PC
84: Input Parameters:
85: + pc - the preconditioner context
86: - levels - number of levels of fill
88: Options Database Key:
89: . -pc_icc_levels <levels> - Sets fill level
91: Level: intermediate
93: Concepts: ICC^setting levels of fill
95: @*/
96: int PCICCSetLevels(PC pc,int levels)
97: {
98: int ierr,(*f)(PC,int);
102: if (levels < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"negative levels");
103: PetscObjectQueryFunction((PetscObject)pc,"PCICCSetLevels_C",(void (**)(void))&f);
104: if (f) {
105: (*f)(pc,levels);
106: }
107: return(0);
108: }
110: /*@
111: PCICCSetFill - Indicate the amount of fill you expect in the factored matrix,
112: where fill = number nonzeros in factor/number nonzeros in original matrix.
114: Collective on PC
116: Input Parameters:
117: + pc - the preconditioner context
118: - fill - amount of expected fill
120: Options Database Key:
121: $ -pc_icc_fill <fill>
123: Note:
124: For sparse matrix factorizations it is difficult to predict how much
125: fill to expect. By running with the option -log_info PETSc will print the
126: actual amount of fill used; allowing you to set the value accurately for
127: future runs. But default PETSc uses a value of 1.0
129: Level: intermediate
131: .keywords: PC, set, factorization, direct, fill
133: .seealso: PCLUSetFill()
134: @*/
135: int PCICCSetFill(PC pc,PetscReal fill)
136: {
137: int ierr,(*f)(PC,PetscReal);
141: if (fill < 1.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Fill factor cannot be less than 1.0");
142: PetscObjectQueryFunction((PetscObject)pc,"PCICCSetFill_C",(void (**)(void))&f);
143: if (f) {
144: (*f)(pc,fill);
145: }
146: return(0);
147: }
149: static int PCSetup_ICC(PC pc)
150: {
151: PC_ICC *icc = (PC_ICC*)pc->data;
152: IS perm,cperm;
153: int ierr;
156: MatGetOrdering(pc->pmat,icc->ordering,&perm,&cperm);
158: if (!pc->setupcalled) {
159: MatICCFactorSymbolic(pc->pmat,perm,icc->fill,icc->levels,&icc->fact);
160: } else if (pc->flag != SAME_NONZERO_PATTERN) {
161: MatDestroy(icc->fact);
162: MatICCFactorSymbolic(pc->pmat,perm,icc->fill,icc->levels,&icc->fact);
163: }
164: ISDestroy(cperm);
165: ISDestroy(perm);
166: MatCholeskyFactorNumeric(pc->pmat,&icc->fact);
167: return(0);
168: }
170: static int PCDestroy_ICC(PC pc)
171: {
172: PC_ICC *icc = (PC_ICC*)pc->data;
173: int ierr;
176: if (icc->fact) {MatDestroy(icc->fact);}
177: PetscStrfree(icc->ordering);
178: PetscFree(icc);
179: return(0);
180: }
182: static int PCApply_ICC(PC pc,Vec x,Vec y)
183: {
184: PC_ICC *icc = (PC_ICC*)pc->data;
185: int ierr;
188: MatSolve(icc->fact,x,y);
189: return(0);
190: }
192: static int PCApplySymmetricLeft_ICC(PC pc,Vec x,Vec y)
193: {
194: int ierr;
195: PC_ICC *icc = (PC_ICC*)pc->data;
198: MatForwardSolve(icc->fact,x,y);
199: return(0);
200: }
202: static int PCApplySymmetricRight_ICC(PC pc,Vec x,Vec y)
203: {
204: int ierr;
205: PC_ICC *icc = (PC_ICC*)pc->data;
208: MatBackwardSolve(icc->fact,x,y);
209: return(0);
210: }
212: static int PCGetFactoredMatrix_ICC(PC pc,Mat *mat)
213: {
214: PC_ICC *icc = (PC_ICC*)pc->data;
217: *mat = icc->fact;
218: return(0);
219: }
221: static int PCSetFromOptions_ICC(PC pc)
222: {
223: PC_ICC *icc = (PC_ICC*)pc->data;
224: char tname[256];
225: PetscTruth flg;
226: int ierr;
227: PetscFList ordlist;
230: MatOrderingRegisterAll(PETSC_NULL);
231: PetscOptionsHead("ICC Options");
232: PetscOptionsInt("-pc_icc_levels","levels of fill","PCICCSetLevels",icc->levels,&icc->levels,&flg);
233: PetscOptionsReal("-pc_icc_fill","Expected fill in factorization","PCICCSetFill",icc->fill,&icc->fill,&flg);
234: MatGetOrderingList(&ordlist);
235: PetscOptionsList("-pc_icc_mat_ordering_type","Reorder to reduce nonzeros in ICC","PCICCSetMatOrdering",ordlist,icc->ordering,tname,256,&flg);
236: if (flg) {
237: PCICCSetMatOrdering(pc,tname);
238: }
239: PetscOptionsTail();
240: return(0);
241: }
243: EXTERN_C_BEGIN
244: int PCCreate_ICC(PC pc)
245: {
246: int ierr;
247: PC_ICC *icc;
250: PetscNew(PC_ICC,&icc);
251: PetscLogObjectMemory(pc,sizeof(PC_ICC));
253: icc->fact = 0;
254: PetscStrallocpy(MATORDERING_NATURAL,&icc->ordering);
255: icc->levels = 0;
256: icc->fill = 1.0;
257: icc->implctx = 0;
258: pc->data = (void*)icc;
260: pc->ops->apply = PCApply_ICC;
261: pc->ops->setup = PCSetup_ICC;
262: pc->ops->destroy = PCDestroy_ICC;
263: pc->ops->setfromoptions = PCSetFromOptions_ICC;
264: pc->ops->view = 0;
265: pc->ops->getfactoredmatrix = PCGetFactoredMatrix_ICC;
266: pc->ops->applysymmetricleft = PCApplySymmetricLeft_ICC;
267: pc->ops->applysymmetricright = PCApplySymmetricRight_ICC;
269: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCICCSetLevels_C","PCICCSetLevels_ICC",
270: PCICCSetLevels_ICC);
271: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCICCSetFill_C","PCICCSetFill_ICC",
272: PCICCSetFill_ICC);
273: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCICCSetMatOrdering_C","PCICCSetMatOrdering_ICC",
274: PCICCSetMatOrdering_ICC);
275: return(0);
276: }
277: EXTERN_C_END