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: #undef __FUNCT__  
 12: int PCICCSetMatOrdering_ICC(PC pc,MatOrderingType ordering)
 13: {
 14:   PC_ICC *dir = (PC_ICC*)pc->data;
 15:   int    ierr;
 16: 
 18:   PetscStrfree(dir->ordering);
 19:   PetscStrallocpy(ordering,&dir->ordering);
 20:   return(0);
 21: }
 22: EXTERN_C_END

 24: EXTERN_C_BEGIN
 25: #undef __FUNCT__  
 27: int PCICCSetFill_ICC(PC pc,PetscReal fill)
 28: {
 29:   PC_ICC *dir;

 32:   dir       = (PC_ICC*)pc->data;
 33:   dir->fill = fill;
 34:   return(0);
 35: }
 36: EXTERN_C_END

 38: EXTERN_C_BEGIN
 39: #undef __FUNCT__  
 41: int PCICCSetLevels_ICC(PC pc,int levels)
 42: {
 43:   PC_ICC *icc;

 46:   icc = (PC_ICC*)pc->data;
 47:   icc->levels = levels;
 48:   return(0);
 49: }
 50: EXTERN_C_END

 52: #undef __FUNCT__  
 54: /*@
 55:     PCICCSetMatOrdering - Sets the ordering routine (to reduce fill) to 
 56:     be used it the ICC factorization.

 58:     Collective on PC

 60:     Input Parameters:
 61: +   pc - the preconditioner context
 62: -   ordering - the matrix ordering name, for example, MATORDERING_ND or MATORDERING_RCM

 64:     Options Database Key:
 65: .   -pc_icc_mat_ordering_type <nd,rcm,...> - Sets ordering routine

 67:     Level: intermediate

 69: .seealso: PCLUSetMatOrdering()

 71: .keywords: PC, ICC, set, matrix, reordering

 73: @*/
 74: int PCICCSetMatOrdering(PC pc,MatOrderingType ordering)
 75: {
 76:   int ierr,(*f)(PC,MatOrderingType);

 80:   PetscObjectQueryFunction((PetscObject)pc,"PCICCSetMatOrdering_C",(void (**)(void))&f);
 81:   if (f) {
 82:     (*f)(pc,ordering);
 83:   }
 84:   return(0);
 85: }

 87: #undef __FUNCT__  
 89: /*@
 90:    PCICCSetLevels - Sets the number of levels of fill to use.

 92:    Collective on PC

 94:    Input Parameters:
 95: +  pc - the preconditioner context
 96: -  levels - number of levels of fill

 98:    Options Database Key:
 99: .  -pc_icc_levels <levels> - Sets fill level

101:    Level: intermediate

103:    Concepts: ICC^setting levels of fill

105: @*/
106: int PCICCSetLevels(PC pc,int levels)
107: {
108:   int ierr,(*f)(PC,int);

112:   if (levels < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"negative levels");
113:   PetscObjectQueryFunction((PetscObject)pc,"PCICCSetLevels_C",(void (**)(void))&f);
114:   if (f) {
115:     (*f)(pc,levels);
116:   }
117:   return(0);
118: }

120: #undef __FUNCT__  
122: /*@
123:    PCICCSetFill - Indicate the amount of fill you expect in the factored matrix,
124:    where fill = number nonzeros in factor/number nonzeros in original matrix.

126:    Collective on PC

128:    Input Parameters:
129: +  pc - the preconditioner context
130: -  fill - amount of expected fill

132:    Options Database Key:
133: $  -pc_icc_fill <fill>

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

141:    Level: intermediate

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

145: .seealso: PCLUSetFill()
146: @*/
147: int PCICCSetFill(PC pc,PetscReal fill)
148: {
149:   int ierr,(*f)(PC,PetscReal);

153:   if (fill < 1.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Fill factor cannot be less than 1.0");
154:   PetscObjectQueryFunction((PetscObject)pc,"PCICCSetFill_C",(void (**)(void))&f);
155:   if (f) {
156:     (*f)(pc,fill);
157:   }
158:   return(0);
159: }

161: #undef __FUNCT__  
163: static int PCSetup_ICC(PC pc)
164: {
165:   PC_ICC *icc = (PC_ICC*)pc->data;
166:   IS     perm,cperm;
167:   int    ierr;

170:   MatGetOrdering(pc->pmat,icc->ordering,&perm,&cperm);

172:   if (!pc->setupcalled) {
173:     MatICCFactorSymbolic(pc->pmat,perm,icc->fill,icc->levels,&icc->fact);
174:   } else if (pc->flag != SAME_NONZERO_PATTERN) {
175:     MatDestroy(icc->fact);
176:     MatICCFactorSymbolic(pc->pmat,perm,icc->fill,icc->levels,&icc->fact);
177:   }
178:   ISDestroy(cperm);
179:   ISDestroy(perm);
180:   MatCholeskyFactorNumeric(pc->pmat,&icc->fact);
181:   return(0);
182: }

184: #undef __FUNCT__  
186: static int PCDestroy_ICC(PC pc)
187: {
188:   PC_ICC *icc = (PC_ICC*)pc->data;
189:   int    ierr;

192:   if (icc->fact) {MatDestroy(icc->fact);}
193:   PetscStrfree(icc->ordering);
194:   PetscFree(icc);
195:   return(0);
196: }

198: #undef __FUNCT__  
200: static int PCApply_ICC(PC pc,Vec x,Vec y)
201: {
202:   PC_ICC *icc = (PC_ICC*)pc->data;
203:   int    ierr;

206:   MatSolve(icc->fact,x,y);
207:   return(0);
208: }

210: #undef __FUNCT__  
212: static int PCApplySymmetricLeft_ICC(PC pc,Vec x,Vec y)
213: {
214:   int    ierr;
215:   PC_ICC *icc = (PC_ICC*)pc->data;

218:   MatForwardSolve(icc->fact,x,y);
219:   return(0);
220: }

222: #undef __FUNCT__  
224: static int PCApplySymmetricRight_ICC(PC pc,Vec x,Vec y)
225: {
226:   int    ierr;
227:   PC_ICC *icc = (PC_ICC*)pc->data;

230:   MatBackwardSolve(icc->fact,x,y);
231:   return(0);
232: }

234: #undef __FUNCT__  
236: static int PCGetFactoredMatrix_ICC(PC pc,Mat *mat)
237: {
238:   PC_ICC *icc = (PC_ICC*)pc->data;

241:   *mat = icc->fact;
242:   return(0);
243: }

245: #undef __FUNCT__  
247: static int PCSetFromOptions_ICC(PC pc)
248: {
249:   PC_ICC     *icc = (PC_ICC*)pc->data;
250:   char       tname[256];
251:   PetscTruth flg;
252:   int        ierr;
253:   PetscFList ordlist;

256:   MatOrderingRegisterAll(PETSC_NULL);
257:   PetscOptionsHead("ICC Options");
258:     PetscOptionsInt("-pc_icc_levels","levels of fill","PCICCSetLevels",icc->levels,&icc->levels,&flg);
259:     PetscOptionsReal("-pc_icc_fill","Expected fill in factorization","PCICCSetFill",icc->fill,&icc->fill,&flg);
260:     MatGetOrderingList(&ordlist);
261:     PetscOptionsList("-pc_icc_mat_ordering_type","Reorder to reduce nonzeros in ICC","PCICCSetMatOrdering",ordlist,icc->ordering,tname,256,&flg);
262:     if (flg) {
263:       PCICCSetMatOrdering(pc,tname);
264:     }
265:   PetscOptionsTail();
266:   return(0);
267: }

269: #undef __FUNCT__  
271: static int PCView_ICC(PC pc,PetscViewer viewer)
272: {
273:   PC_ICC     *icc = (PC_ICC*)pc->data;
274:   int        ierr;
275:   PetscTruth isstring,isascii;

278:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
279:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
280:   if (isascii) {
281:     if (icc->levels == 1) {
282:         PetscViewerASCIIPrintf(viewer,"  ICC: %d level of filln",icc->levels);
283:     } else {
284:         PetscViewerASCIIPrintf(viewer,"  ICC: %d levels of filln",icc->levels);
285:     }
286:     PetscViewerASCIIPrintf(viewer,"  ICC: max fill ratio allocated %gn",icc->fill);
287:   } else if (isstring) {
288:     PetscViewerStringSPrintf(viewer," lvls=%d",icc->levels);
289:   } else {
290:     SETERRQ1(1,"Viewer type %s not supported for PCICC",((PetscObject)viewer)->type_name);
291:   }
292:   return(0);
293: }

295: EXTERN_C_BEGIN
296: #undef __FUNCT__  
298: int PCCreate_ICC(PC pc)
299: {
300:   int    ierr;
301:   PC_ICC *icc;

304:   PetscNew(PC_ICC,&icc);
305:   PetscLogObjectMemory(pc,sizeof(PC_ICC));

307:   icc->fact                  = 0;
308:   PetscStrallocpy(MATORDERING_NATURAL,&icc->ordering);
309:   icc->levels                  = 0;
310:   icc->fill               = 1.0;
311:   icc->implctx            = 0;
312:   pc->data                  = (void*)icc;

314:   pc->ops->apply               = PCApply_ICC;
315:   pc->ops->setup               = PCSetup_ICC;
316:   pc->ops->destroy               = PCDestroy_ICC;
317:   pc->ops->setfromoptions      = PCSetFromOptions_ICC;
318:   pc->ops->view                = PCView_ICC;
319:   pc->ops->getfactoredmatrix   = PCGetFactoredMatrix_ICC;
320:   pc->ops->applysymmetricleft  = PCApplySymmetricLeft_ICC;
321:   pc->ops->applysymmetricright = PCApplySymmetricRight_ICC;

323:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCICCSetLevels_C","PCICCSetLevels_ICC",
324:                     PCICCSetLevels_ICC);
325:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCICCSetFill_C","PCICCSetFill_ICC",
326:                     PCICCSetFill_ICC);
327:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCICCSetMatOrdering_C","PCICCSetMatOrdering_ICC",
328:                     PCICCSetMatOrdering_ICC);
329:   return(0);
330: }
331: EXTERN_C_END