Actual source code: hypre.c

  1: #define PETSCKSP_DLL

  3: /* 
  4:    Provides an interface to the LLNL package hypre
  5: */

  7: /* Must use hypre 2.0.0 or more recent. */

 9:  #include private/pcimpl.h
 11: #include "HYPRE.h"
 12: #include "HYPRE_parcsr_ls.h"
 13: #include "_hypre_parcsr_mv.h"
 14: #include "_hypre_IJ_mv.h"

 17: EXTERN PetscErrorCode MatHYPRE_IJMatrixCreate(Mat,HYPRE_IJMatrix*);
 18: EXTERN PetscErrorCode MatHYPRE_IJMatrixCopy(Mat,HYPRE_IJMatrix);
 19: EXTERN PetscErrorCode MatHYPRE_IJMatrixFastCopy(Mat,HYPRE_IJMatrix);
 20: EXTERN PetscErrorCode VecHYPRE_IJVectorCreate(Vec,HYPRE_IJVector*);

 22: /* 
 23:    Private context (data structure) for the  preconditioner.  
 24: */
 25: typedef struct {
 26:   HYPRE_Solver       hsolver;
 27:   HYPRE_IJMatrix     ij;
 28:   HYPRE_IJVector     b,x;

 30:   PetscErrorCode     (*destroy)(HYPRE_Solver);
 31:   PetscErrorCode     (*solve)(HYPRE_Solver,HYPRE_ParCSRMatrix,HYPRE_ParVector,HYPRE_ParVector);
 32:   PetscErrorCode     (*setup)(HYPRE_Solver,HYPRE_ParCSRMatrix,HYPRE_ParVector,HYPRE_ParVector);
 33: 
 34:   MPI_Comm           comm_hypre;
 35:   char              *hypre_type;

 37:   /* options for Pilut and BoomerAMG*/
 38:   int                maxiter;
 39:   double             tol;

 41:   /* options for Pilut */
 42:   int                factorrowsize;

 44:   /* options for ParaSails */
 45:   int                nlevels;
 46:   double             threshhold;
 47:   double             filter;
 48:   int                sym;
 49:   double             loadbal;
 50:   int                logging;
 51:   int                ruse;
 52:   int                symt;

 54:   /* options for Euclid */
 55:   PetscTruth         bjilu;
 56:   int                levels;

 58:   /* options for Euclid and BoomerAMG */
 59:   PetscTruth         printstatistics;

 61:   /* options for BoomerAMG */
 62:   int                cycletype;
 63:   int                maxlevels;
 64:   double             strongthreshold;
 65:   double             maxrowsum;
 66:   int                gridsweeps[3];
 67:   int                coarsentype;
 68:   int                measuretype;
 69:   int                relaxtype[3];
 70:   double             relaxweight;
 71:   double             outerrelaxweight;
 72:   int                relaxorder;
 73:   double             truncfactor;
 74:   PetscTruth         applyrichardson;
 75:   int                pmax;
 76:   int                interptype;
 77:   int                agg_nl;
 78:   int                agg_num_paths;
 79:   int                nodal_coarsen;
 80:   PetscTruth         nodal_relax;
 81:   int                nodal_relax_levels;
 82: } PC_HYPRE;


 87: static PetscErrorCode PCSetUp_HYPRE(PC pc)
 88: {
 89:   PC_HYPRE           *jac = (PC_HYPRE*)pc->data;
 90:   PetscErrorCode     ierr;
 91:   HYPRE_ParCSRMatrix hmat;
 92:   HYPRE_ParVector    bv,xv;
 93:   PetscInt           bs;
 94:   int                hierr;

 97:   if (!jac->hypre_type) {
 98:     PCHYPRESetType(pc,"boomeramg");
 99:   }

101:   if (pc->setupcalled) {
102:     /* always destroy the old matrix and create a new memory;
103:        hope this does not churn the memory too much. The problem
104:        is I do not know if it is possible to put the matrix back to
105:        its initial state so that we can directly copy the values 
106:        the second time through. */
107:     HYPRE_IJMatrixDestroy(jac->ij);
108:     jac->ij = 0;
109:   }

111:   if (!jac->ij) { /* create the matrix the first time through */
112:     MatHYPRE_IJMatrixCreate(pc->pmat,&jac->ij);
113:   }
114:   if (!jac->b) { /* create the vectors the first time through */
115:     Vec x,b;
116:     MatGetVecs(pc->pmat,&x,&b);
117:     VecHYPRE_IJVectorCreate(x,&jac->x);
118:     VecHYPRE_IJVectorCreate(b,&jac->b);
119:     VecDestroy(x);
120:     VecDestroy(b);
121:   }

123:   /* special case for BoomerAMG */
124:   if (jac->setup == HYPRE_BoomerAMGSetup) {
125:     MatGetBlockSize(pc->pmat,&bs);
126:     if (bs > 1) {
127:       HYPRE_BoomerAMGSetNumFunctions(jac->hsolver,bs);
128:     }
129:   };
130:   MatHYPRE_IJMatrixCopy(pc->pmat,jac->ij);
131:   HYPRE_IJMatrixGetObject(jac->ij,(void**)&hmat);
132:   HYPRE_IJVectorGetObject(jac->b,(void**)&bv);
133:   HYPRE_IJVectorGetObject(jac->x,(void**)&xv);
134:   h(*jac->setup)(jac->hsolver,hmat,bv,xv);
135:   if (hierr) SETERRQ1(PETSC_ERR_LIB,"Error in HYPRE setup, error code %d",hierr);
136:   return(0);
137: }

139: /*
140:     Replaces the address where the HYPRE vector points to its data with the address of
141:   PETSc's data. Saves the old address so it can be reset when we are finished with it.
142:   Allows use to get the data into a HYPRE vector without the cost of memcopies 
143: */
144: #define HYPREReplacePointer(b,newvalue,savedvalue) {\
145:    hypre_ParVector *par_vector   = (hypre_ParVector *)hypre_IJVectorObject(((hypre_IJVector*)b));\
146:    hypre_Vector    *local_vector = hypre_ParVectorLocalVector(par_vector);\
147:    savedvalue         = local_vector->data;\
148:    local_vector->data = newvalue;}

152: static PetscErrorCode PCApply_HYPRE(PC pc,Vec b,Vec x)
153: {
154:   PC_HYPRE           *jac = (PC_HYPRE*)pc->data;
155:   PetscErrorCode     ierr;
156:   HYPRE_ParCSRMatrix hmat;
157:   PetscScalar        *bv,*xv;
158:   HYPRE_ParVector    jbv,jxv;
159:   PetscScalar        *sbv,*sxv;
160:   int                hierr;

163:   if (!jac->applyrichardson) {VecSet(x,0.0);}
164:   VecGetArray(b,&bv);
165:   VecGetArray(x,&xv);
166:   HYPREReplacePointer(jac->b,bv,sbv);
167:   HYPREReplacePointer(jac->x,xv,sxv);

169:   HYPRE_IJMatrixGetObject(jac->ij,(void**)&hmat);
170:   HYPRE_IJVectorGetObject(jac->b,(void**)&jbv);
171:   HYPRE_IJVectorGetObject(jac->x,(void**)&jxv);
172:   h(*jac->solve)(jac->hsolver,hmat,jbv,jxv);

174:   /*if (hierr && (hierr != HYPRE_ERROR_CONV || jac->solve != HYPRE_BoomerAMGSolve))SETERRQ1(PETSC_ERR_LIB,"Error in HYPRE solver, error code %d",hierr);
175:    */
176:  /* error code of HYPRE_ERROR_CONV means convergence not achieved - if
177:     the tolerance is set to 0.0 (the default), a convergence error will
178:     not occur (so we may not want to overide the conv. error here?*/
179:  if (hierr && hierr != HYPRE_ERROR_CONV)
180:   {
181:      SETERRQ1(PETSC_ERR_LIB,"Error in HYPRE solver, error code %d",hierr);
182:   }
183:  if (hierr) hypre__global_error = 0;
184: 

186:   HYPREReplacePointer(jac->b,sbv,bv);
187:   HYPREReplacePointer(jac->x,sxv,xv);
188:   VecRestoreArray(x,&xv);
189:   VecRestoreArray(b,&bv);
190:   return(0);
191: }

195: static PetscErrorCode PCDestroy_HYPRE(PC pc)
196: {
197:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;

201:   if (jac->ij) { HYPRE_IJMatrixDestroy(jac->ij); }
202:   if (jac->b)  { HYPRE_IJVectorDestroy(jac->b);  }
203:   if (jac->x)  { HYPRE_IJVectorDestroy(jac->x);  }
204:   if (jac->destroy) { (*jac->destroy)(jac->hsolver); }
205:   PetscStrfree(jac->hypre_type);
206:   if (jac->comm_hypre != MPI_COMM_NULL) { MPI_Comm_free(&(jac->comm_hypre));}
207:   PetscFree(jac);

209:   PetscObjectChangeTypeName((PetscObject)pc,0);
210:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCHYPRESetType_C","",PETSC_NULL);
211:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCHYPREGetType_C","",PETSC_NULL);
212:   return(0);
213: }

215: /* --------------------------------------------------------------------------------------------*/
218: static PetscErrorCode PCSetFromOptions_HYPRE_Pilut(PC pc)
219: {
220:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
222:   PetscTruth     flag;

225:   PetscOptionsHead("HYPRE Pilut Options");
226:   PetscOptionsInt("-pc_hypre_pilut_maxiter","Number of iterations","None",jac->maxiter,&jac->maxiter,&flag);
227:   if (flag) {
228:     HYPRE_ParCSRPilutSetMaxIter(jac->hsolver,jac->maxiter);
229:   }
230:   PetscOptionsReal("-pc_hypre_pilut_tol","Drop tolerance","None",jac->tol,&jac->tol,&flag);
231:   if (flag) {
232:     HYPRE_ParCSRPilutSetDropTolerance(jac->hsolver,jac->tol);
233:   }
234:   PetscOptionsInt("-pc_hypre_pilut_factorrowsize","FactorRowSize","None",jac->factorrowsize,&jac->factorrowsize,&flag);
235:   if (flag) {
236:     HYPRE_ParCSRPilutSetFactorRowSize(jac->hsolver,jac->factorrowsize);
237:   }
238:   PetscOptionsTail();
239:   return(0);
240: }

244: static PetscErrorCode PCView_HYPRE_Pilut(PC pc,PetscViewer viewer)
245: {
246:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
248:   PetscTruth     iascii;

251:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
252:   if (iascii) {
253:     PetscViewerASCIIPrintf(viewer,"  HYPRE Pilut preconditioning\n");
254:     if (jac->maxiter != PETSC_DEFAULT) {
255:       PetscViewerASCIIPrintf(viewer,"  HYPRE Pilut: maximum number of iterations %d\n",jac->maxiter);
256:     } else {
257:       PetscViewerASCIIPrintf(viewer,"  HYPRE Pilut: default maximum number of iterations \n");
258:     }
259:     if (jac->tol != PETSC_DEFAULT) {
260:       PetscViewerASCIIPrintf(viewer,"  HYPRE Pilut: drop tolerance %G\n",jac->tol);
261:     } else {
262:       PetscViewerASCIIPrintf(viewer,"  HYPRE Pilut: default drop tolerance \n");
263:     }
264:     if (jac->factorrowsize != PETSC_DEFAULT) {
265:       PetscViewerASCIIPrintf(viewer,"  HYPRE Pilut: factor row size %d\n",jac->factorrowsize);
266:     } else {
267:       PetscViewerASCIIPrintf(viewer,"  HYPRE Pilut: default factor row size \n");
268:     }
269:   }
270:   return(0);
271: }

273: /* --------------------------------------------------------------------------------------------*/
276: static PetscErrorCode PCSetFromOptions_HYPRE_Euclid(PC pc)
277: {
278:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
280:   PetscTruth     flag;
281:   char           *args[8],levels[16];
282:   PetscInt       cnt = 0;

285:   PetscOptionsHead("HYPRE Euclid Options");
286:   PetscOptionsInt("-pc_hypre_euclid_levels","Number of levels of fill ILU(k)","None",jac->levels,&jac->levels,&flag);
287:   if (flag) {
288:     if (jac->levels < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Number of levels %d must be nonegative",jac->levels);
289:     sprintf(levels,"%d",jac->levels);
290:     args[cnt++] = (char*)"-level"; args[cnt++] = levels;
291:   }
292:   PetscOptionsTruth("-pc_hypre_euclid_bj","Use block Jacobi ILU(k)","None",jac->bjilu,&jac->bjilu,PETSC_NULL);
293:   if (jac->bjilu) {
294:     args[cnt++] =(char*) "-bj"; args[cnt++] = (char*)"1";
295:   }
296: 
297:   PetscOptionsTruth("-pc_hypre_euclid_print_statistics","Print statistics","None",jac->printstatistics,&jac->printstatistics,PETSC_NULL);
298:   if (jac->printstatistics) {
299:     args[cnt++] = (char*)"-eu_stats"; args[cnt++] = (char*)"1";
300:     args[cnt++] = (char*)"-eu_mem"; args[cnt++] = (char*)"1";
301:   }
302:   PetscOptionsTail();
303:   if (cnt) {
304:     HYPRE_EuclidSetParams(jac->hsolver,cnt,args);
305:   }
306:   return(0);
307: }

311: static PetscErrorCode PCView_HYPRE_Euclid(PC pc,PetscViewer viewer)
312: {
313:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
315:   PetscTruth     iascii;

318:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
319:   if (iascii) {
320:     PetscViewerASCIIPrintf(viewer,"  HYPRE Euclid preconditioning\n");
321:     PetscViewerASCIIPrintf(viewer,"  HYPRE Euclid: number of levels %d\n",jac->levels);
322:     if (jac->bjilu) {
323:       PetscViewerASCIIPrintf(viewer,"  HYPRE Euclid: Using block Jacobi ILU instead of parallel ILU\n");
324:     }
325:   }
326:   return(0);
327: }

329: /* --------------------------------------------------------------------------------------------*/

333: static PetscErrorCode PCApplyTranspose_HYPRE_BoomerAMG(PC pc,Vec b,Vec x)
334: {
335:   PC_HYPRE           *jac = (PC_HYPRE*)pc->data;
336:   PetscErrorCode     ierr;
337:   HYPRE_ParCSRMatrix hmat;
338:   PetscScalar        *bv,*xv;
339:   HYPRE_ParVector    jbv,jxv;
340:   PetscScalar        *sbv,*sxv;
341:   int                hierr;

344:   VecSet(x,0.0);
345:   VecGetArray(b,&bv);
346:   VecGetArray(x,&xv);
347:   HYPREReplacePointer(jac->b,bv,sbv);
348:   HYPREReplacePointer(jac->x,xv,sxv);

350:   HYPRE_IJMatrixGetObject(jac->ij,(void**)&hmat);
351:   HYPRE_IJVectorGetObject(jac->b,(void**)&jbv);
352:   HYPRE_IJVectorGetObject(jac->x,(void**)&jxv);
353: 
354:   hHYPRE_BoomerAMGSolveT(jac->hsolver,hmat,jbv,jxv);
355:   /* error code of 1 in BoomerAMG merely means convergence not achieved */
356:   if (hierr && (hierr != 1)) SETERRQ1(PETSC_ERR_LIB,"Error in HYPRE solver, error code %d",hierr);
357:   if (hierr) hypre__global_error = 0;
358: 
359:   HYPREReplacePointer(jac->b,sbv,bv);
360:   HYPREReplacePointer(jac->x,sxv,xv);
361:   VecRestoreArray(x,&xv);
362:   VecRestoreArray(b,&bv);
363:   return(0);
364: }

366: static const char *HYPREBoomerAMGCycleType[]   = {"","V","W"};
367: static const char *HYPREBoomerAMGCoarsenType[] = {"CLJP","Ruge-Stueben","","modifiedRuge-Stueben","","","Falgout", "", "PMIS", "", "HMIS"};
368: static const char *HYPREBoomerAMGMeasureType[] = {"local","global"};
369: static const char *HYPREBoomerAMGRelaxType[]   = {"Jacobi","sequential-Gauss-Seidel","","SOR/Jacobi","backward-SOR/Jacobi","","symmetric-SOR/Jacobi",
370:                                                   "","","Gaussian-elimination"};
371: static const char *HYPREBoomerAMGInterpType[]  = {"classical", "", "", "direct", "multipass", "multipass-wts", "ext+i",
372:                                                   "ext+i-cc", "standard", "standard-wts", "", "", "FF", "FF1"};
375: static PetscErrorCode PCSetFromOptions_HYPRE_BoomerAMG(PC pc)
376: {
377:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
379:   int            n,indx;
380:   PetscTruth     flg, tmp_truth;
381:   double         tmpdbl, twodbl[2];

384:   PetscOptionsHead("HYPRE BoomerAMG Options");
385:   PetscOptionsEList("-pc_hypre_boomeramg_cycle_type","Cycle type","None",HYPREBoomerAMGCycleType,2,HYPREBoomerAMGCycleType[jac->cycletype],&indx,&flg);
386:   if (flg) {
387:     jac->cycletype = indx;
388:     HYPRE_BoomerAMGSetCycleType(jac->hsolver,jac->cycletype);
389:   }
390:   PetscOptionsInt("-pc_hypre_boomeramg_max_levels","Number of levels (of grids) allowed","None",jac->maxlevels,&jac->maxlevels,&flg);
391:   if (flg) {
392:     if (jac->maxlevels < 2) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Number of levels %d must be at least two",jac->maxlevels);
393:     HYPRE_BoomerAMGSetMaxLevels(jac->hsolver,jac->maxlevels);
394:   }
395:   PetscOptionsInt("-pc_hypre_boomeramg_max_iter","Maximum iterations used PER hypre call","None",jac->maxiter,&jac->maxiter,&flg);
396:   if (flg) {
397:     if (jac->maxiter < 1) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Number of iterations %d must be at least one",jac->maxiter);
398:     HYPRE_BoomerAMGSetMaxIter(jac->hsolver,jac->maxiter);
399:   }
400:   PetscOptionsScalar("-pc_hypre_boomeramg_tol","Convergence tolerance PER hypre call (0.0 = use a fixed number of iterations)","None",jac->tol,&jac->tol,&flg);
401:   if (flg) {
402:     if (jac->tol < 0.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Tolerance %G must be greater than or equal to zero",jac->tol);
403:     HYPRE_BoomerAMGSetTol(jac->hsolver,jac->tol);
404:   }

406:   PetscOptionsScalar("-pc_hypre_boomeramg_truncfactor","Truncation factor for interpolation (0=no truncation)","None",jac->truncfactor,&jac->truncfactor,&flg);
407:   if (flg) {
408:     if (jac->truncfactor < 0.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Truncation factor %G must be great than or equal zero",jac->truncfactor);
409:     HYPRE_BoomerAMGSetTruncFactor(jac->hsolver,jac->truncfactor);
410:   }

412:  PetscOptionsInt("-pc_hypre_boomeramg_P_max","Max elements per row for interpolation operator ( 0=unlimited )","None",jac->pmax,&jac->pmax,&flg);
413:   if (flg) {
414:     if (jac->pmax < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"P_max %G must be greater than or equal to zero",jac->pmax);
415:     HYPRE_BoomerAMGSetPMaxElmts(jac->hsolver,jac->pmax);
416:   }

418:  PetscOptionsInt("-pc_hypre_boomeramg_agg_nl","Number of levels of aggressive coarsening","None",jac->agg_nl,&jac->agg_nl,&flg);
419:   if (flg) {
420:      if (jac->agg_nl < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Number of levels %G must be greater than or equal to zero",jac->agg_nl);

422:      HYPRE_BoomerAMGSetAggNumLevels(jac->hsolver,jac->agg_nl);
423:   }


426:  PetscOptionsInt("-pc_hypre_boomeramg_agg_num_paths","Number of paths for aggressive coarsening","None",jac->agg_num_paths,&jac->agg_num_paths,&flg);
427:   if (flg) {
428:      if (jac->agg_num_paths < 1) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Number of paths %G must be greater than or equal to 1",jac->agg_num_paths);

430:      HYPRE_BoomerAMGSetNumPaths(jac->hsolver,jac->agg_num_paths);
431:   }


434:   PetscOptionsScalar("-pc_hypre_boomeramg_strong_threshold","Threshold for being strongly connected","None",jac->strongthreshold,&jac->strongthreshold,&flg);
435:   if (flg) {
436:     if (jac->strongthreshold < 0.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Strong threshold %G must be great than or equal zero",jac->strongthreshold);
437:     HYPRE_BoomerAMGSetStrongThreshold(jac->hsolver,jac->strongthreshold);
438:   }
439:   PetscOptionsScalar("-pc_hypre_boomeramg_max_row_sum","Maximum row sum","None",jac->maxrowsum,&jac->maxrowsum,&flg);
440:   if (flg) {
441:     if (jac->maxrowsum < 0.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Maximum row sum %G must be greater than zero",jac->maxrowsum);
442:     if (jac->maxrowsum > 1.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Maximum row sum %G must be less than or equal one",jac->maxrowsum);
443:     HYPRE_BoomerAMGSetMaxRowSum(jac->hsolver,jac->maxrowsum);
444:   }

446:   /* Grid sweeps */
447:   PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_all","Number of sweeps for the up and down grid levels","None",jac->gridsweeps[0],&indx,&flg);
448:   if (flg) {
449:     HYPRE_BoomerAMGSetNumSweeps(jac->hsolver,indx);
450:     /* modify the jac structure so we can view the updated options with PC_View */
451:     jac->gridsweeps[0] = indx;
452:     jac->gridsweeps[1] = indx;
453:     /*defaults coarse to 1 */
454:     jac->gridsweeps[2] = 1;
455:   }

457:   PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_down","Number of sweeps for the down cycles","None",jac->gridsweeps[0], &indx ,&flg);
458:   if (flg) {
459:     HYPRE_BoomerAMGSetCycleNumSweeps(jac->hsolver,indx, 1);
460:     jac->gridsweeps[0] = indx;
461:   }
462:   PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_up","Number of sweeps for the up cycles","None",jac->gridsweeps[1],&indx,&flg);
463:   if (flg) {
464:     HYPRE_BoomerAMGSetCycleNumSweeps(jac->hsolver,indx, 2);
465:     jac->gridsweeps[1] = indx;
466:   }
467:   PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_coarse","Number of sweeps for the coarse level","None",jac->gridsweeps[2],&indx,&flg);
468:   if (flg) {
469:     HYPRE_BoomerAMGSetCycleNumSweeps(jac->hsolver,indx, 3);
470:     jac->gridsweeps[2] = indx;
471:   }

473:   /* Relax type */
474:   PetscOptionsEList("-pc_hypre_boomeramg_relax_type_all","Relax type for the up and down cycles","None",HYPREBoomerAMGRelaxType,10,HYPREBoomerAMGRelaxType[6],&indx,&flg);
475:   if (flg) {
476:     jac->relaxtype[0] = jac->relaxtype[1]  = indx;
477:     HYPRE_BoomerAMGSetRelaxType(jac->hsolver, indx);
478:     /* by default, coarse type set to 9 */
479:     jac->relaxtype[2] = 9;
480: 
481:   }
482:   PetscOptionsEList("-pc_hypre_boomeramg_relax_type_down","Relax type for the down cycles","None",HYPREBoomerAMGRelaxType,10,HYPREBoomerAMGRelaxType[6],&indx,&flg);
483:   if (flg) {
484:     jac->relaxtype[0] = indx;
485:     HYPRE_BoomerAMGSetCycleRelaxType(jac->hsolver, indx, 1);
486:   }
487:   PetscOptionsEList("-pc_hypre_boomeramg_relax_type_up","Relax type for the up cycles","None",HYPREBoomerAMGRelaxType,10,HYPREBoomerAMGRelaxType[6],&indx,&flg);
488:   if (flg) {
489:     jac->relaxtype[1] = indx;
490:     HYPRE_BoomerAMGSetCycleRelaxType(jac->hsolver, indx, 2);
491:   }
492:   PetscOptionsEList("-pc_hypre_boomeramg_relax_type_coarse","Relax type on coarse grid","None",HYPREBoomerAMGRelaxType,10,HYPREBoomerAMGRelaxType[9],&indx,&flg);
493:   if (flg) {
494:     jac->relaxtype[2] = indx;
495:     HYPRE_BoomerAMGSetCycleRelaxType(jac->hsolver, indx, 3);
496:   }

498:   /* Relaxation Weight */
499:   PetscOptionsReal("-pc_hypre_boomeramg_relax_weight_all","Relaxation weight for all levels (0 = hypre estimates, -k = determined with k CG steps)","None",jac->relaxweight, &tmpdbl ,&flg);
500:   if (flg) {
501:     HYPRE_BoomerAMGSetRelaxWt(jac->hsolver,tmpdbl);
502:     jac->relaxweight = tmpdbl;
503:   }

505:   n=2;
506:   twodbl[0] = twodbl[1] = 1.0;
507:   PetscOptionsRealArray("-pc_hypre_boomeramg_relax_weight_level","Set the relaxation weight for a particular level (weight,level)","None",twodbl, &n, &flg);
508:   if (flg) {
509:     if (n == 2) {
510:       indx =  (int)PetscAbsReal(twodbl[1]);
511:       HYPRE_BoomerAMGSetLevelRelaxWt(jac->hsolver,twodbl[0],indx);
512:     } else {
513:       SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Relax weight level: you must provide 2 values separated by a comma (and no space), you provided %d",n);
514:     }
515:   }

517:   /* Outer relaxation Weight */
518:   PetscOptionsReal("-pc_hypre_boomeramg_outer_relax_weight_all","Outer relaxation weight for all levels ( -k = determined with k CG steps)","None",jac->outerrelaxweight, &tmpdbl ,&flg);
519:   if (flg) {
520:     HYPRE_BoomerAMGSetOuterWt( jac->hsolver, tmpdbl);
521:     jac->outerrelaxweight = tmpdbl;
522:   }

524:   n=2;
525:   twodbl[0] = twodbl[1] = 1.0;
526:   PetscOptionsRealArray("-pc_hypre_boomeramg_outer_relax_weight_level","Set the outer relaxation weight for a particular level (weight,level)","None",twodbl, &n, &flg);
527:   if (flg) {
528:     if (n == 2) {
529:       indx =  (int)PetscAbsReal(twodbl[1]);
530:       HYPRE_BoomerAMGSetLevelOuterWt( jac->hsolver, twodbl[0], indx);
531:     } else {
532:       SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Relax weight outer level: You must provide 2 values separated by a comma (and no space), you provided %d",n);
533:     }
534:   }

536:   /* the Relax Order */
537:   /* PetscOptionsName("-pc_hypre_boomeramg_no_CF", "Do not use CF-relaxation", "None", &flg); */
538:   PetscOptionsTruth( "-pc_hypre_boomeramg_no_CF", "Do not use CF-relaxation", "None", PETSC_FALSE, &tmp_truth, &flg);

540:   if (flg) {
541:     jac->relaxorder = 0;
542:     HYPRE_BoomerAMGSetRelaxOrder(jac->hsolver, jac->relaxorder);
543:   }
544:   PetscOptionsEList("-pc_hypre_boomeramg_measure_type","Measure type","None",HYPREBoomerAMGMeasureType,2,HYPREBoomerAMGMeasureType[0],&indx,&flg);
545:   if (flg) {
546:     jac->measuretype = indx;
547:     HYPRE_BoomerAMGSetMeasureType(jac->hsolver,jac->measuretype);
548:   }
549:   /* update list length 3/07 */
550:   PetscOptionsEList("-pc_hypre_boomeramg_coarsen_type","Coarsen type","None",HYPREBoomerAMGCoarsenType,11,HYPREBoomerAMGCoarsenType[6],&indx,&flg);
551:   if (flg) {
552:     jac->coarsentype = indx;
553:     HYPRE_BoomerAMGSetCoarsenType(jac->hsolver,jac->coarsentype);
554:   }
555: 
556:   /* new 3/07 */
557:   PetscOptionsEList("-pc_hypre_boomeramg_interp_type","Interpolation type","None",HYPREBoomerAMGInterpType,14,HYPREBoomerAMGInterpType[0],&indx,&flg);
558:   if (flg) {
559:     jac->interptype = indx;
560:     HYPRE_BoomerAMGSetInterpType(jac->hsolver,jac->interptype);
561:   }


564:   PetscOptionsName("-pc_hypre_boomeramg_print_statistics","Print statistics","None",&flg);
565:   if (flg) {
566:     int level=3;
567:     jac->printstatistics = PETSC_TRUE;
568:     PetscOptionsInt("-pc_hypre_boomeramg_print_statistics","Print statistics","None",level,&level,PETSC_NULL);
569:     HYPRE_BoomerAMGSetPrintLevel(jac->hsolver,level);
570:   }

572:   PetscOptionsName("-pc_hypre_boomeramg_print_debug","Print debug information","None",&flg);
573:   if (flg) {
574:     int level=3;
575:     jac->printstatistics = PETSC_TRUE;
576:     PetscOptionsInt("-pc_hypre_boomeramg_print_debug","Print debug information","None",level,&level,PETSC_NULL);
577:     HYPRE_BoomerAMGSetDebugFlag(jac->hsolver,level);
578:   }

580:   PetscOptionsTruth( "-pc_hypre_boomeramg_nodal_coarsen", "HYPRE_BoomerAMGSetNodal()", "None", PETSC_FALSE, &tmp_truth, &flg);
581:   if (flg && tmp_truth) {
582:     jac->nodal_coarsen = 1;
583:     HYPRE_BoomerAMGSetNodal(jac->hsolver,1);
584:   }

586:   PetscOptionsTruth( "-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None", PETSC_FALSE, &tmp_truth, &flg);
587:   if (flg && tmp_truth) {
588:     PetscInt tmp_int;
589:     PetscOptionsInt( "-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None",jac->nodal_relax_levels,&tmp_int,&flg);
590:     if (flg) jac->nodal_relax_levels = tmp_int;
591:     HYPRE_BoomerAMGSetSmoothType(jac->hsolver,6);
592:     HYPRE_BoomerAMGSetDomainType(jac->hsolver,1);
593:     HYPRE_BoomerAMGSetOverlap(jac->hsolver,0);
594:     HYPRE_BoomerAMGSetSmoothNumLevels(jac->hsolver,jac->nodal_relax_levels);
595:   }

597:   PetscOptionsTail();
598:   return(0);
599: }

603: static PetscErrorCode PCApplyRichardson_HYPRE_BoomerAMG(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its,PetscInt *outits,PCRichardsonConvergedReason *reason)
604: {
605:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
607:   int            oits;

610:   HYPRE_BoomerAMGSetMaxIter(jac->hsolver,its*jac->maxiter);
611:   HYPRE_BoomerAMGSetTol(jac->hsolver,rtol);
612:   jac->applyrichardson = PETSC_TRUE;
613:   PCApply_HYPRE(pc,b,y);
614:   jac->applyrichardson = PETSC_FALSE;
615:   HYPRE_BoomerAMGGetNumIterations(jac->hsolver,&oits);
616:   *outits = oits;
617:   if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
618:   else             *reason = PCRICHARDSON_CONVERGED_RTOL;
619:   HYPRE_BoomerAMGSetTol(jac->hsolver,jac->tol);
620:   HYPRE_BoomerAMGSetMaxIter(jac->hsolver,jac->maxiter);
621:   return(0);
622: }


627: static PetscErrorCode PCView_HYPRE_BoomerAMG(PC pc,PetscViewer viewer)
628: {
629:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
631:   PetscTruth     iascii;

634:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
635:   if (iascii) {
636:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG preconditioning\n");
637:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Cycle type %s\n",HYPREBoomerAMGCycleType[jac->cycletype]);
638:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Maximum number of levels %d\n",jac->maxlevels);
639:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Maximum number of iterations PER hypre call %d\n",jac->maxiter);
640:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Convergence tolerance PER hypre call %G\n",jac->tol);
641:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Threshold for strong coupling %G\n",jac->strongthreshold);
642:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Interpolation truncation factor %G\n",jac->truncfactor);
643:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Interpolation: max elements per row %d\n",jac->pmax);
644:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Number of levels of aggressive coarsening %d\n",jac->agg_nl);
645:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Number of paths for aggressive coarsening %d\n",jac->agg_num_paths);
646: 
647:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Maximum row sums %G\n",jac->maxrowsum);

649:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Sweeps down         %d\n",jac->gridsweeps[0]);
650:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Sweeps up           %d\n",jac->gridsweeps[1]);
651:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Sweeps on coarse    %d\n",jac->gridsweeps[2]);

653:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Relax down          %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[0]]);
654:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Relax up            %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[1]]);
655:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Relax on coarse     %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[2]]);

657:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Relax weight  (all)      %G\n",jac->relaxweight);
658:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Outer relax weight (all) %G\n",jac->outerrelaxweight);

660:     if (jac->relaxorder) {
661:       PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Using CF-relaxation\n");
662:     } else {
663:       PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Not using CF-relaxation\n");
664:     }
665:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Measure type        %s\n",HYPREBoomerAMGMeasureType[jac->measuretype]);
666:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Coarsen type        %s\n",HYPREBoomerAMGCoarsenType[jac->coarsentype]);
667:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG: Interpolation type  %s\n",HYPREBoomerAMGInterpType[jac->interptype]);
668:     if (jac->nodal_coarsen) {
669:       PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Using nodal coarsening (with HYPRE_BOOMERAMGSetNodal())\n");
670:     }
671:     if (jac->nodal_relax) {
672:       PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Using nodal relaxation via Schwarz smoothing on levels %d\n",jac->nodal_relax_levels);
673:     }
674:   }
675:   return(0);
676: }

678: /* --------------------------------------------------------------------------------------------*/
681: static PetscErrorCode PCSetFromOptions_HYPRE_ParaSails(PC pc)
682: {
683:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
685:   int            indx;
686:   PetscTruth     flag;
687:   const char     *symtlist[] = {"nonsymmetric","SPD","nonsymmetric,SPD"};

690:   PetscOptionsHead("HYPRE ParaSails Options");
691:   PetscOptionsInt("-pc_hypre_parasails_nlevels","Number of number of levels","None",jac->nlevels,&jac->nlevels,0);
692:   PetscOptionsReal("-pc_hypre_parasails_thresh","Threshold","None",jac->threshhold,&jac->threshhold,&flag);
693:   if (flag) {
694:     HYPRE_ParaSailsSetParams(jac->hsolver,jac->threshhold,jac->nlevels);
695:   }

697:   PetscOptionsReal("-pc_hypre_parasails_filter","filter","None",jac->filter,&jac->filter,&flag);
698:   if (flag) {
699:     HYPRE_ParaSailsSetFilter(jac->hsolver,jac->filter);
700:   }

702:   PetscOptionsReal("-pc_hypre_parasails_loadbal","Load balance","None",jac->loadbal,&jac->loadbal,&flag);
703:   if (flag) {
704:     HYPRE_ParaSailsSetLoadbal(jac->hsolver,jac->loadbal);
705:   }

707:   PetscOptionsTruth("-pc_hypre_parasails_logging","Print info to screen","None",(PetscTruth)jac->logging,(PetscTruth*)&jac->logging,&flag);
708:   if (flag) {
709:     HYPRE_ParaSailsSetLogging(jac->hsolver,jac->logging);
710:   }

712:   PetscOptionsTruth("-pc_hypre_parasails_reuse","Reuse nonzero pattern in preconditioner","None",(PetscTruth)jac->ruse,(PetscTruth*)&jac->ruse,&flag);
713:   if (flag) {
714:     HYPRE_ParaSailsSetReuse(jac->hsolver,jac->ruse);
715:   }

717:   PetscOptionsEList("-pc_hypre_parasails_sym","Symmetry of matrix and preconditioner","None",symtlist,3,symtlist[0],&indx,&flag);
718:   if (flag) {
719:     jac->symt = indx;
720:     HYPRE_ParaSailsSetSym(jac->hsolver,jac->symt);
721:   }

723:   PetscOptionsTail();
724:   return(0);
725: }

729: static PetscErrorCode PCView_HYPRE_ParaSails(PC pc,PetscViewer viewer)
730: {
731:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
733:   PetscTruth     iascii;
734:   const char     *symt = 0;;

737:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
738:   if (iascii) {
739:     PetscViewerASCIIPrintf(viewer,"  HYPRE ParaSails preconditioning\n");
740:     PetscViewerASCIIPrintf(viewer,"  HYPRE ParaSails: nlevels %d\n",jac->nlevels);
741:     PetscViewerASCIIPrintf(viewer,"  HYPRE ParaSails: threshold %G\n",jac->threshhold);
742:     PetscViewerASCIIPrintf(viewer,"  HYPRE ParaSails: filter %G\n",jac->filter);
743:     PetscViewerASCIIPrintf(viewer,"  HYPRE ParaSails: load balance %G\n",jac->loadbal);
744:     PetscViewerASCIIPrintf(viewer,"  HYPRE ParaSails: reuse nonzero structure %s\n",PetscTruths[jac->ruse]);
745:     PetscViewerASCIIPrintf(viewer,"  HYPRE ParaSails: print info to screen %s\n",PetscTruths[jac->logging]);
746:     if (!jac->symt) {
747:       symt = "nonsymmetric matrix and preconditioner";
748:     } else if (jac->symt == 1) {
749:       symt = "SPD matrix and preconditioner";
750:     } else if (jac->symt == 2) {
751:       symt = "nonsymmetric matrix but SPD preconditioner";
752:     } else {
753:       SETERRQ1(PETSC_ERR_ARG_WRONG,"Unknown HYPRE ParaSails symmetric option %d",jac->symt);
754:     }
755:     PetscViewerASCIIPrintf(viewer,"  HYPRE ParaSails: %s\n",symt);
756:   }
757:   return(0);
758: }
759: /* ---------------------------------------------------------------------------------*/

764: PetscErrorCode  PCHYPREGetType_HYPRE(PC pc,const char *name[])
765: {
766:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;

769:   *name = jac->hypre_type;
770:   return(0);
771: }

777: PetscErrorCode  PCHYPRESetType_HYPRE(PC pc,const char name[])
778: {
779:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
781:   PetscTruth     flag;

784:   if (jac->hypre_type) {
785:     PetscStrcmp(jac->hypre_type,name,&flag);
786:     if (!flag) SETERRQ(PETSC_ERR_ORDER,"Cannot reset the HYPRE preconditioner type once it has been set");
787:     return(0);
788:   } else {
789:     PetscStrallocpy(name, &jac->hypre_type);
790:   }
791: 
792:   jac->maxiter            = PETSC_DEFAULT;
793:   jac->tol                = PETSC_DEFAULT;
794:   jac->printstatistics    = PetscLogPrintInfo;

796:   PetscStrcmp("pilut",jac->hypre_type,&flag);
797:   if (flag) {
798:     HYPRE_ParCSRPilutCreate(jac->comm_hypre,&jac->hsolver);
799:     pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Pilut;
800:     pc->ops->view           = PCView_HYPRE_Pilut;
801:     jac->destroy            = HYPRE_ParCSRPilutDestroy;
802:     jac->setup              = HYPRE_ParCSRPilutSetup;
803:     jac->solve              = HYPRE_ParCSRPilutSolve;
804:     jac->factorrowsize      = PETSC_DEFAULT;
805:     return(0);
806:   }
807:   PetscStrcmp("parasails",jac->hypre_type,&flag);
808:   if (flag) {
809:     HYPRE_ParaSailsCreate(jac->comm_hypre,&jac->hsolver);
810:     pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ParaSails;
811:     pc->ops->view           = PCView_HYPRE_ParaSails;
812:     jac->destroy            = HYPRE_ParaSailsDestroy;
813:     jac->setup              = HYPRE_ParaSailsSetup;
814:     jac->solve              = HYPRE_ParaSailsSolve;
815:     /* initialize */
816:     jac->nlevels     = 1;
817:     jac->threshhold  = .1;
818:     jac->filter      = .1;
819:     jac->loadbal     = 0;
820:     if (PetscLogPrintInfo) {
821:       jac->logging     = (int) PETSC_TRUE;
822:     } else {
823:       jac->logging     = (int) PETSC_FALSE;
824:     }
825:     jac->ruse = (int) PETSC_FALSE;
826:     jac->symt   = 0;
827:     HYPRE_ParaSailsSetParams(jac->hsolver,jac->threshhold,jac->nlevels);
828:     HYPRE_ParaSailsSetFilter(jac->hsolver,jac->filter);
829:     HYPRE_ParaSailsSetLoadbal(jac->hsolver,jac->loadbal);
830:     HYPRE_ParaSailsSetLogging(jac->hsolver,jac->logging);
831:     HYPRE_ParaSailsSetReuse(jac->hsolver,jac->ruse);
832:     HYPRE_ParaSailsSetSym(jac->hsolver,jac->symt);
833:     return(0);
834:   }
835:   PetscStrcmp("euclid",jac->hypre_type,&flag);
836:   if (flag) {
837:     HYPRE_EuclidCreate(jac->comm_hypre,&jac->hsolver);
838:     pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Euclid;
839:     pc->ops->view           = PCView_HYPRE_Euclid;
840:     jac->destroy            = HYPRE_EuclidDestroy;
841:     jac->setup              = HYPRE_EuclidSetup;
842:     jac->solve              = HYPRE_EuclidSolve;
843:     /* initialization */
844:     jac->bjilu              = PETSC_FALSE;
845:     jac->levels             = 1;
846:     return(0);
847:   }
848:   PetscStrcmp("boomeramg",jac->hypre_type,&flag);
849:   if (flag) {
850:     HYPRE_BoomerAMGCreate(&jac->hsolver);
851:     pc->ops->setfromoptions  = PCSetFromOptions_HYPRE_BoomerAMG;
852:     pc->ops->view            = PCView_HYPRE_BoomerAMG;
853:     pc->ops->applytranspose  = PCApplyTranspose_HYPRE_BoomerAMG;
854:     pc->ops->applyrichardson = PCApplyRichardson_HYPRE_BoomerAMG;
855:     jac->destroy             = HYPRE_BoomerAMGDestroy;
856:     jac->setup               = HYPRE_BoomerAMGSetup;
857:     jac->solve               = HYPRE_BoomerAMGSolve;
858:     jac->applyrichardson     = PETSC_FALSE;
859:     /* these defaults match the hypre defaults */
860:     jac->cycletype        = 1;
861:     jac->maxlevels        = 25;
862:     jac->maxiter          = 1;
863:     jac->tol              = 0.0; /* tolerance of zero indicates use as preconditioner (suppresses convergence errors) */
864:     jac->truncfactor      = 0.0;
865:     jac->strongthreshold  = .25;
866:     jac->maxrowsum        = .9;
867:     jac->coarsentype      = 6;
868:     jac->measuretype      = 0;
869:     jac->gridsweeps[0]    = jac->gridsweeps[1] = jac->gridsweeps[2] = 1;
870:     jac->relaxtype[0]     = jac->relaxtype[1] = 6; /* Defaults to SYMMETRIC since in PETSc we are using a a PC - most likely with CG */
871:     jac->relaxtype[2]     = 9; /*G.E. */
872:     jac->relaxweight      = 1.0;
873:     jac->outerrelaxweight = 1.0;
874:     jac->relaxorder       = 1;
875:     jac->interptype       = 0;
876:     jac->agg_nl           = 0;
877:     jac->pmax             = 0;
878:     jac->truncfactor      = 0.0;
879:     jac->agg_num_paths    = 1;

881:     jac->nodal_coarsen    = 0;
882:     jac->nodal_relax      = PETSC_FALSE;
883:     jac->nodal_relax_levels = 1;
884:     HYPRE_BoomerAMGSetCycleType(jac->hsolver,jac->cycletype);
885:     HYPRE_BoomerAMGSetMaxLevels(jac->hsolver,jac->maxlevels);
886:     HYPRE_BoomerAMGSetMaxIter(jac->hsolver,jac->maxiter);
887:     HYPRE_BoomerAMGSetTol(jac->hsolver,jac->tol);
888:     HYPRE_BoomerAMGSetTruncFactor(jac->hsolver,jac->truncfactor);
889:     HYPRE_BoomerAMGSetStrongThreshold(jac->hsolver,jac->strongthreshold);
890:     HYPRE_BoomerAMGSetMaxRowSum(jac->hsolver,jac->maxrowsum);
891:     HYPRE_BoomerAMGSetCoarsenType(jac->hsolver,jac->coarsentype);
892:     HYPRE_BoomerAMGSetMeasureType(jac->hsolver,jac->measuretype);
893:     HYPRE_BoomerAMGSetRelaxOrder(jac->hsolver, jac->relaxorder);
894:     HYPRE_BoomerAMGSetInterpType(jac->hsolver,jac->interptype);
895:     HYPRE_BoomerAMGSetAggNumLevels(jac->hsolver,jac->agg_nl);
896:     HYPRE_BoomerAMGSetPMaxElmts(jac->hsolver,jac->pmax);
897:     HYPRE_BoomerAMGSetNumPaths(jac->hsolver,jac->agg_num_paths);
898:     HYPRE_BoomerAMGSetRelaxType(jac->hsolver, jac->relaxtype[0]);  /*defaults coarse to 9*/
899:     HYPRE_BoomerAMGSetNumSweeps(jac->hsolver, jac->gridsweeps[0]); /*defaults coarse to 1 */
900:     return(0);
901:   }
902:   PetscStrfree(jac->hypre_type);
903:   jac->hypre_type = PETSC_NULL;
904:   SETERRQ1(PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown HYPRE preconditioner %s; Choices are pilut, parasails, euclid, boomeramg",name);
905:   return(0);
906: }

909: /*
910:     It only gets here if the HYPRE type has not been set before the call to 
911:    ...SetFromOptions() which actually is most of the time
912: */
915: static PetscErrorCode PCSetFromOptions_HYPRE(PC pc)
916: {
918:   int            indx;
919:   const char     *type[] = {"pilut","parasails","boomeramg","euclid"};
920:   PetscTruth     flg;

923:   PetscOptionsHead("HYPRE preconditioner options");
924:   PetscOptionsEList("-pc_hypre_type","HYPRE preconditioner type","PCHYPRESetType",type,4,"boomeramg",&indx,&flg);
925:   if (flg) {
926:     PCHYPRESetType_HYPRE(pc,type[indx]);
927:   } else {
928:     PCHYPRESetType_HYPRE(pc,"boomeramg");
929:   }
930:   if (pc->ops->setfromoptions) {
931:     pc->ops->setfromoptions(pc);
932:   }
933:   PetscOptionsTail();
934:   return(0);
935: }

939: /*@C
940:      PCHYPRESetType - Sets which hypre preconditioner you wish to use

942:    Input Parameters:
943: +     pc - the preconditioner context
944: -     name - either  pilut, parasails, boomeramg, euclid

946:    Options Database Keys:
947:    -pc_hypre_type - One of pilut, parasails, boomeramg, euclid
948:  
949:    Level: intermediate

951: .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PC,
952:            PCHYPRE

954: @*/
955: PetscErrorCode  PCHYPRESetType(PC pc,const char name[])
956: {
957:   PetscErrorCode ierr,(*f)(PC,const char[]);

962:   PetscObjectQueryFunction((PetscObject)pc,"PCHYPRESetType_C",(void (**)(void))&f);
963:   if (f) {
964:     (*f)(pc,name);
965:   }
966:   return(0);
967: }

971: /*@C
972:      PCHYPREGetType - Gets which hypre preconditioner you are using

974:    Input Parameter:
975: .     pc - the preconditioner context

977:    Output Parameter:
978: .     name - either  pilut, parasails, boomeramg, euclid

980:    Level: intermediate

982: .seealso:  PCCreate(), PCHYPRESetType(), PCType (for list of available types), PC,
983:            PCHYPRE

985: @*/
986: PetscErrorCode  PCHYPREGetType(PC pc,const char *name[])
987: {
988:   PetscErrorCode ierr,(*f)(PC,const char*[]);

993:   PetscObjectQueryFunction((PetscObject)pc,"PCHYPREGetType_C",(void (**)(void))&f);
994:   if (f) {
995:     (*f)(pc,name);
996:   }
997:   return(0);
998: }

1000: /*MC
1001:      PCHYPRE - Allows you to use the matrix element based preconditioners in the LLNL package hypre

1003:    Options Database Keys:
1004: +   -pc_hypre_type - One of pilut, parasails, boomeramg, euclid
1005: -   Too many others to list, run with -pc_type hypre -pc_hypre_type XXX -help to see options for the XXX
1006:           preconditioner
1007:  
1008:    Level: intermediate

1010:    Notes: Apart from pc_hypre_type (for which there is PCHYPRESetType()),
1011:           the many hypre options can ONLY be set via the options database (e.g. the command line
1012:           or with PetscOptionsSetValue(), there are no functions to set them)

1014:           The options -pc_hypre_boomeramg_max_iter and -pc_hypre_boomeramg_rtol refer to the number of iterations
1015:           (V-cycles) and tolerance that boomeramg does EACH time it is called. So for example, if 
1016:           -pc_hypre_boomeramg_max_iter is set to 2 then 2-V-cycles are being used to define the preconditioner 
1017:           (-pc_hypre_boomeramg_rtol should be set to 0.0 - the default - to strictly use a fixed number of 
1018:           iterations per hypre call). -ksp_max_it and -ksp_rtol STILL determine the total number of iterations 
1019:           and tolerance for the Krylov solver. For example, if -pc_hypre_boomeramg_max_iter is 2 and -ksp_max_it is 10 
1020:           then AT MOST twenty V-cycles of boomeramg will be called.

1022:            Note that the option -pc_hypre_boomeramg_relax_type_all defaults to symmetric relaxation 
1023:            (symmetric-SOR/Jacobi), which is required for Krylov solvers like CG that expect symmetry.  
1024:            Otherwise, you may want to use -pc_hypre_boomeramg_relax_type_all SOR/Jacobi.
1025:           If you wish to use BoomerAMG WITHOUT a Krylov method use -ksp_type richardson NOT -ksp_type preonly
1026:           and use -ksp_max_it to control the number of V-cycles.
1027:           (see the PETSc FAQ.html at the PETSc website under the Documentation tab).

1029:           2007-02-03 Using HYPRE-1.11.1b, the routine HYPRE_BoomerAMGSolveT and the option
1030:           -pc_hypre_parasails_reuse were failing with SIGSEGV. Dalcin L.

1032: .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PC,
1033:            PCHYPRESetType()

1035: M*/

1040: PetscErrorCode  PCCreate_HYPRE(PC pc)
1041: {
1042:   PC_HYPRE       *jac;

1046:   PetscNewLog(pc,PC_HYPRE,&jac);
1047:   pc->data                 = jac;
1048:   pc->ops->destroy         = PCDestroy_HYPRE;
1049:   pc->ops->setfromoptions  = PCSetFromOptions_HYPRE;
1050:   pc->ops->setup           = PCSetUp_HYPRE;
1051:   pc->ops->apply           = PCApply_HYPRE;
1052:   jac->comm_hypre          = MPI_COMM_NULL;
1053:   jac->hypre_type          = PETSC_NULL;
1054:   /* duplicate communicator for hypre */
1055:   MPI_Comm_dup(((PetscObject)pc)->comm,&(jac->comm_hypre));
1056:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCHYPRESetType_C","PCHYPRESetType_HYPRE",PCHYPRESetType_HYPRE);
1057:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCHYPREGetType_C","PCHYPREGetType_HYPRE",PCHYPREGetType_HYPRE);
1058:   return(0);
1059: }