Actual source code: hypre.c

petsc-3.8.2 2017-11-09
Report Typos and Errors

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

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

  8:  #include <petsc/private/pcimpl.h>
  9: /* this include is needed ONLY to allow access to the private data inside the Mat object specific to hypre */
 10:  #include <petsc/private/matimpl.h>
 11:  #include <../src/vec/vec/impls/hypre/vhyp.h>
 12:  #include <../src/mat/impls/hypre/mhypre.h>
 13:  #include <../src/dm/impls/da/hypre/mhyp.h>
 14: #include <_hypre_parcsr_ls.h>

 16: static PetscBool cite = PETSC_FALSE;
 17: static const char hypreCitation[] = "@manual{hypre-web-page,\n  title  = {{\\sl hypre}: High Performance Preconditioners},\n  organization = {Lawrence Livermore National Laboratory},\n  note  = {\\url{http://www.llnl.gov/CASC/hypre/}}\n}\n";

 19: /*
 20:    Private context (data structure) for the  preconditioner.
 21: */
 22: typedef struct {
 23:   HYPRE_Solver   hsolver;
 24:   Mat            hpmat; /* MatHYPRE */

 26:   HYPRE_Int (*destroy)(HYPRE_Solver);
 27:   HYPRE_Int (*solve)(HYPRE_Solver,HYPRE_ParCSRMatrix,HYPRE_ParVector,HYPRE_ParVector);
 28:   HYPRE_Int (*setup)(HYPRE_Solver,HYPRE_ParCSRMatrix,HYPRE_ParVector,HYPRE_ParVector);

 30:   MPI_Comm comm_hypre;
 31:   char     *hypre_type;

 33:   /* options for Pilut and BoomerAMG*/
 34:   PetscInt maxiter;
 35:   double   tol;

 37:   /* options for Pilut */
 38:   PetscInt factorrowsize;

 40:   /* options for ParaSails */
 41:   PetscInt nlevels;
 42:   double   threshhold;
 43:   double   filter;
 44:   PetscInt sym;
 45:   double   loadbal;
 46:   PetscInt logging;
 47:   PetscInt ruse;
 48:   PetscInt symt;

 50:   /* options for BoomerAMG */
 51:   PetscBool printstatistics;

 53:   /* options for BoomerAMG */
 54:   PetscInt  cycletype;
 55:   PetscInt  maxlevels;
 56:   double    strongthreshold;
 57:   double    maxrowsum;
 58:   PetscInt  gridsweeps[3];
 59:   PetscInt  coarsentype;
 60:   PetscInt  measuretype;
 61:   PetscInt  smoothtype;
 62:   PetscInt  smoothnumlevels;
 63:   PetscInt  eu_level;   /* Number of levels for ILU(k) in Euclid */
 64:   double    eu_droptolerance; /* Drop tolerance for ILU(k) in Euclid */
 65:   PetscInt  eu_bj;      /* Defines use of Block Jacobi ILU in Euclid */
 66:   PetscInt  relaxtype[3];
 67:   double    relaxweight;
 68:   double    outerrelaxweight;
 69:   PetscInt  relaxorder;
 70:   double    truncfactor;
 71:   PetscBool applyrichardson;
 72:   PetscInt  pmax;
 73:   PetscInt  interptype;
 74:   PetscInt  agg_nl;
 75:   PetscInt  agg_num_paths;
 76:   PetscInt  nodal_coarsen;
 77:   PetscBool nodal_relax;
 78:   PetscInt  nodal_relax_levels;

 80:   PetscInt  nodal_coarsening;
 81:   PetscInt  vec_interp_variant;
 82:   HYPRE_IJVector  *hmnull;
 83:   HYPRE_ParVector *phmnull;  /* near null space passed to hypre */
 84:   PetscInt        n_hmnull;
 85:   Vec             hmnull_constant;
 86:   PetscScalar     **hmnull_hypre_data_array;   /* this is the space in hmnull that was allocated by hypre, it is restored to hypre just before freeing the phmnull vectors */

 88:   /* options for AS (Auxiliary Space preconditioners) */
 89:   PetscInt  as_print;
 90:   PetscInt  as_max_iter;
 91:   PetscReal as_tol;
 92:   PetscInt  as_relax_type;
 93:   PetscInt  as_relax_times;
 94:   PetscReal as_relax_weight;
 95:   PetscReal as_omega;
 96:   PetscInt  as_amg_alpha_opts[5]; /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for vector Poisson (AMS) or Curl problem (ADS) */
 97:   PetscReal as_amg_alpha_theta;   /* AMG strength for vector Poisson (AMS) or Curl problem (ADS) */
 98:   PetscInt  as_amg_beta_opts[5];  /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for scalar Poisson (AMS) or vector Poisson (ADS) */
 99:   PetscReal as_amg_beta_theta;    /* AMG strength for scalar Poisson (AMS) or vector Poisson (ADS)  */
100:   PetscInt  ams_cycle_type;
101:   PetscInt  ads_cycle_type;

103:   /* additional data */
104:   Mat G;             /* MatHYPRE */
105:   Mat C;             /* MatHYPRE */
106:   Mat alpha_Poisson; /* MatHYPRE */
107:   Mat beta_Poisson;  /* MatHYPRE */

109:   /* extra information for AMS */
110:   PetscInt       dim; /* geometrical dimension */
111:   HYPRE_IJVector coords[3];
112:   HYPRE_IJVector constants[3];
113:   Mat            RT_PiFull, RT_Pi[3];
114:   Mat            ND_PiFull, ND_Pi[3];
115:   PetscBool      ams_beta_is_zero;
116:   PetscBool      ams_beta_is_zero_part;
117:   PetscInt       ams_proj_freq;
118: } PC_HYPRE;

120: PetscErrorCode PCHYPREGetSolver(PC pc,HYPRE_Solver *hsolver)
121: {
122:   PC_HYPRE *jac = (PC_HYPRE*)pc->data;

125:   *hsolver = jac->hsolver;
126:   return(0);
127: }

129: static PetscErrorCode PCSetUp_HYPRE(PC pc)
130: {
131:   PC_HYPRE           *jac = (PC_HYPRE*)pc->data;
132:   Mat_HYPRE          *hjac;
133:   HYPRE_ParCSRMatrix hmat;
134:   HYPRE_ParVector    bv,xv;
135:   PetscBool          ishypre;
136:   PetscErrorCode     ierr;

139:   if (!jac->hypre_type) {
140:     PCHYPRESetType(pc,"boomeramg");
141:   }

143:   PetscObjectTypeCompare((PetscObject)pc->pmat,MATHYPRE,&ishypre);
144:   if (!ishypre) {
145:     MatDestroy(&jac->hpmat);
146:     MatConvert(pc->pmat,MATHYPRE,MAT_INITIAL_MATRIX,&jac->hpmat);
147:   } else {
148:     PetscObjectReference((PetscObject)pc->pmat);
149:     MatDestroy(&jac->hpmat);
150:     jac->hpmat = pc->pmat;
151:   }
152:   hjac = (Mat_HYPRE*)(jac->hpmat->data);

154:   /* special case for BoomerAMG */
155:   if (jac->setup == HYPRE_BoomerAMGSetup) {
156:     MatNullSpace    mnull;
157:     PetscBool       has_const;
158:     PetscInt        bs,nvec,i;
159:     const Vec       *vecs;
160:     PetscScalar     *petscvecarray;

162:     MatGetBlockSize(pc->pmat,&bs);
163:     if (bs > 1) PetscStackCallStandard(HYPRE_BoomerAMGSetNumFunctions,(jac->hsolver,bs));
164:     MatGetNearNullSpace(pc->mat, &mnull);
165:     if (mnull) {
166:       MatNullSpaceGetVecs(mnull, &has_const, &nvec, &vecs);
167:       PetscMalloc1(nvec+1,&jac->hmnull);
168:       PetscMalloc1(nvec+1,&jac->hmnull_hypre_data_array);
169:       PetscMalloc1(nvec+1,&jac->phmnull);
170:       for (i=0; i<nvec; i++) {
171:         VecHYPRE_IJVectorCreate(vecs[i],&jac->hmnull[i]);
172:         VecGetArrayRead(vecs[i],(const PetscScalar **)&petscvecarray);
173:         VecHYPRE_ParVectorReplacePointer(jac->hmnull[i],petscvecarray,jac->hmnull_hypre_data_array[i]);
174:         VecRestoreArrayRead(vecs[i],(const PetscScalar **)&petscvecarray);
175:         PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->hmnull[i],(void**)&jac->phmnull[i]));
176:       }
177:       if (has_const) {
178:         MatCreateVecs(pc->pmat,&jac->hmnull_constant,NULL);
179:         VecSet(jac->hmnull_constant,1);
180:         VecNormalize(jac->hmnull_constant,NULL);
181:         VecHYPRE_IJVectorCreate(jac->hmnull_constant,&jac->hmnull[nvec]);
182:         VecGetArrayRead(jac->hmnull_constant,(const PetscScalar **)&petscvecarray);
183:         VecHYPRE_ParVectorReplacePointer(jac->hmnull[nvec],petscvecarray,jac->hmnull_hypre_data_array[nvec]);
184:         VecRestoreArrayRead(jac->hmnull_constant,(const PetscScalar **)&petscvecarray);
185:         PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->hmnull[nvec],(void**)&jac->phmnull[nvec]));
186:         nvec++;
187:       }
188:       PetscStackCallStandard(HYPRE_BoomerAMGSetInterpVectors,(jac->hsolver,nvec,jac->phmnull));
189:       jac->n_hmnull = nvec;
190:     }
191:   }

193:   /* special case for AMS */
194:   if (jac->setup == HYPRE_AMSSetup) {
195:     Mat_HYPRE          *hm;
196:     HYPRE_ParCSRMatrix parcsr;
197:     if (!jac->coords[0] && !jac->constants[0] && !(jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1]))) {
198:       SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE AMS preconditioner needs either the coordinate vectors via PCSetCoordinates() or the edge constant vectors via PCHYPRESetEdgeConstantVectors() or the interpolation matrix via PCHYPRESetInterpolations");
199:     }
200:     if (jac->dim) {
201:       PetscStackCallStandard(HYPRE_AMSSetDimension,(jac->hsolver,jac->dim));
202:     }
203:     if (jac->constants[0]) {
204:       HYPRE_ParVector ozz,zoz,zzo = NULL;
205:       PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->constants[0],(void**)(&ozz)));
206:       PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->constants[1],(void**)(&zoz)));
207:       if (jac->constants[2]) {
208:         PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->constants[2],(void**)(&zzo)));
209:       }
210:       PetscStackCallStandard(HYPRE_AMSSetEdgeConstantVectors,(jac->hsolver,ozz,zoz,zzo));
211:     }
212:     if (jac->coords[0]) {
213:       HYPRE_ParVector coords[3];
214:       coords[0] = NULL;
215:       coords[1] = NULL;
216:       coords[2] = NULL;
217:       if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[0],(void**)(&coords[0])));
218:       if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[1],(void**)(&coords[1])));
219:       if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[2],(void**)(&coords[2])));
220:       PetscStackCallStandard(HYPRE_AMSSetCoordinateVectors,(jac->hsolver,coords[0],coords[1],coords[2]));
221:     }
222:     if (!jac->G) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE AMS preconditioner needs the discrete gradient operator via PCHYPRESetDiscreteGradient");
223:     hm = (Mat_HYPRE*)(jac->G->data);
224:     PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&parcsr)));
225:     PetscStackCallStandard(HYPRE_AMSSetDiscreteGradient,(jac->hsolver,parcsr));
226:     if (jac->alpha_Poisson) {
227:       hm = (Mat_HYPRE*)(jac->alpha_Poisson->data);
228:       PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&parcsr)));
229:       PetscStackCallStandard(HYPRE_AMSSetAlphaPoissonMatrix,(jac->hsolver,parcsr));
230:     }
231:     if (jac->ams_beta_is_zero) {
232:       PetscStackCallStandard(HYPRE_AMSSetBetaPoissonMatrix,(jac->hsolver,NULL));
233:     } else if (jac->beta_Poisson) {
234:       hm = (Mat_HYPRE*)(jac->beta_Poisson->data);
235:       PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&parcsr)));
236:       PetscStackCallStandard(HYPRE_AMSSetBetaPoissonMatrix,(jac->hsolver,parcsr));
237:     }
238:     if (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1])) {
239:       PetscInt           i;
240:       HYPRE_ParCSRMatrix nd_parcsrfull, nd_parcsr[3];
241:       if (jac->ND_PiFull) {
242:         hm = (Mat_HYPRE*)(jac->ND_PiFull->data);
243:         PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&nd_parcsrfull)));
244:       } else {
245:         nd_parcsrfull = NULL;
246:       }
247:       for (i=0;i<3;++i) {
248:         if (jac->ND_Pi[i]) {
249:           hm = (Mat_HYPRE*)(jac->ND_Pi[i]->data);
250:           PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&nd_parcsr[i])));
251:         } else {
252:           nd_parcsr[i] = NULL;
253:         }
254:       }
255:       PetscStackCallStandard(HYPRE_AMSSetInterpolations,(jac->hsolver,nd_parcsrfull,nd_parcsr[0],nd_parcsr[1],nd_parcsr[2]));
256:     }
257:   }
258:   /* special case for ADS */
259:   if (jac->setup == HYPRE_ADSSetup) {
260:     Mat_HYPRE          *hm;
261:     HYPRE_ParCSRMatrix parcsr;
262:     if (!jac->coords[0] && !((jac->RT_PiFull || (jac->RT_Pi[0] && jac->RT_Pi[1])) && (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1])))) {
263:       SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs either the coordinate vectors via PCSetCoordinates() or the interpolation matrices via PCHYPRESetInterpolations");
264:     }
265:     else if (!jac->coords[1] || !jac->coords[2]) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner has been designed for three dimensional problems! For two dimensional problems, use HYPRE AMS instead");
266:     if (!jac->G) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs the discrete gradient operator via PCHYPRESetDiscreteGradient");
267:     if (!jac->C) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs the discrete curl operator via PCHYPRESetDiscreteGradient");
268:     if (jac->coords[0]) {
269:       HYPRE_ParVector coords[3];
270:       coords[0] = NULL;
271:       coords[1] = NULL;
272:       coords[2] = NULL;
273:       if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[0],(void**)(&coords[0])));
274:       if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[1],(void**)(&coords[1])));
275:       if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[2],(void**)(&coords[2])));
276:       PetscStackCallStandard(HYPRE_ADSSetCoordinateVectors,(jac->hsolver,coords[0],coords[1],coords[2]));
277:     }
278:     hm = (Mat_HYPRE*)(jac->G->data);
279:     PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&parcsr)));
280:     PetscStackCallStandard(HYPRE_ADSSetDiscreteGradient,(jac->hsolver,parcsr));
281:     hm = (Mat_HYPRE*)(jac->C->data);
282:     PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&parcsr)));
283:     PetscStackCallStandard(HYPRE_ADSSetDiscreteCurl,(jac->hsolver,parcsr));
284:     if ((jac->RT_PiFull || (jac->RT_Pi[0] && jac->RT_Pi[1])) && (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1]))) {
285:       PetscInt           i;
286:       HYPRE_ParCSRMatrix rt_parcsrfull, rt_parcsr[3];
287:       HYPRE_ParCSRMatrix nd_parcsrfull, nd_parcsr[3];
288:       if (jac->RT_PiFull) {
289:         hm = (Mat_HYPRE*)(jac->RT_PiFull->data);
290:         PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&rt_parcsrfull)));
291:       } else {
292:         rt_parcsrfull = NULL;
293:       }
294:       for (i=0;i<3;++i) {
295:         if (jac->RT_Pi[i]) {
296:           hm = (Mat_HYPRE*)(jac->RT_Pi[i]->data);
297:           PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&rt_parcsr[i])));
298:         } else {
299:           rt_parcsr[i] = NULL;
300:         }
301:       }
302:       if (jac->ND_PiFull) {
303:         hm = (Mat_HYPRE*)(jac->ND_PiFull->data);
304:         PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&nd_parcsrfull)));
305:       } else {
306:         nd_parcsrfull = NULL;
307:       }
308:       for (i=0;i<3;++i) {
309:         if (jac->ND_Pi[i]) {
310:           hm = (Mat_HYPRE*)(jac->ND_Pi[i]->data);
311:           PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&nd_parcsr[i])));
312:         } else {
313:           nd_parcsr[i] = NULL;
314:         }
315:       }
316:       PetscStackCallStandard(HYPRE_ADSSetInterpolations,(jac->hsolver,rt_parcsrfull,rt_parcsr[0],rt_parcsr[1],rt_parcsr[2],nd_parcsrfull,nd_parcsr[0],nd_parcsr[1],nd_parcsr[2]));
317:     }
318:   }
319:   PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hjac->ij,(void**)&hmat));
320:   PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->b,(void**)&bv));
321:   PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->x,(void**)&xv));
322:   PetscStackCall("HYPRE_SetupXXX",(*jac->setup)(jac->hsolver,hmat,bv,xv););
323:   return(0);
324: }

326: static PetscErrorCode PCApply_HYPRE(PC pc,Vec b,Vec x)
327: {
328:   PC_HYPRE           *jac = (PC_HYPRE*)pc->data;
329:   Mat_HYPRE          *hjac = (Mat_HYPRE*)(jac->hpmat->data);
330:   PetscErrorCode     ierr;
331:   HYPRE_ParCSRMatrix hmat;
332:   PetscScalar        *xv;
333:   const PetscScalar  *bv,*sbv;
334:   HYPRE_ParVector    jbv,jxv;
335:   PetscScalar        *sxv;
336:   PetscInt           hierr;

339:   PetscCitationsRegister(hypreCitation,&cite);
340:   if (!jac->applyrichardson) {VecSet(x,0.0);}
341:   VecGetArrayRead(b,&bv);
342:   VecGetArray(x,&xv);
343:   VecHYPRE_ParVectorReplacePointer(hjac->b,(PetscScalar*)bv,sbv);
344:   VecHYPRE_ParVectorReplacePointer(hjac->x,xv,sxv);

346:   PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hjac->ij,(void**)&hmat));
347:   PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->b,(void**)&jbv));
348:   PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->x,(void**)&jxv));
349:   PetscStackCall("Hypre solve",h(*jac->solve)(jac->hsolver,hmat,jbv,jxv);
350:   if (hierr && hierr != HYPRE_ERROR_CONV) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"Error in HYPRE solver, error code %d",hierr);
351:   if (hierr) hypre__global_error = 0;);

353:   if (jac->setup == HYPRE_AMSSetup && jac->ams_beta_is_zero_part) {
354:     PetscStackCallStandard(HYPRE_AMSProjectOutGradients,(jac->hsolver,jxv));
355:   }
356:   VecHYPRE_ParVectorReplacePointer(hjac->b,(PetscScalar*)sbv,bv);
357:   VecHYPRE_ParVectorReplacePointer(hjac->x,sxv,xv);
358:   VecRestoreArray(x,&xv);
359:   VecRestoreArrayRead(b,&bv);
360:   return(0);
361: }

363: static PetscErrorCode PCReset_HYPRE(PC pc)
364: {
365:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;

369:   MatDestroy(&jac->hpmat);
370:   MatDestroy(&jac->G);
371:   MatDestroy(&jac->C);
372:   MatDestroy(&jac->alpha_Poisson);
373:   MatDestroy(&jac->beta_Poisson);
374:   MatDestroy(&jac->RT_PiFull);
375:   MatDestroy(&jac->RT_Pi[0]);
376:   MatDestroy(&jac->RT_Pi[1]);
377:   MatDestroy(&jac->RT_Pi[2]);
378:   MatDestroy(&jac->ND_PiFull);
379:   MatDestroy(&jac->ND_Pi[0]);
380:   MatDestroy(&jac->ND_Pi[1]);
381:   MatDestroy(&jac->ND_Pi[2]);
382:   if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[0])); jac->coords[0] = NULL;
383:   if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[1])); jac->coords[1] = NULL;
384:   if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[2])); jac->coords[2] = NULL;
385:   if (jac->constants[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[0])); jac->constants[0] = NULL;
386:   if (jac->constants[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[1])); jac->constants[1] = NULL;
387:   if (jac->constants[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[2])); jac->constants[2] = NULL;
388:   if (jac->n_hmnull && jac->hmnull) {
389:     PetscInt                 i;
390:     PETSC_UNUSED PetscScalar *petscvecarray;

392:     for (i=0; i<jac->n_hmnull; i++) {
393:       VecHYPRE_ParVectorReplacePointer(jac->hmnull[i],jac->hmnull_hypre_data_array[i],petscvecarray);
394:       PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->hmnull[i]));
395:     }
396:     PetscFree(jac->hmnull);
397:     PetscFree(jac->hmnull_hypre_data_array);
398:     PetscFree(jac->phmnull);
399:     VecDestroy(&jac->hmnull_constant);
400:   }
401:   jac->ams_beta_is_zero = PETSC_FALSE;
402:   jac->dim = 0;
403:   return(0);
404: }

406: static PetscErrorCode PCDestroy_HYPRE(PC pc)
407: {
408:   PC_HYPRE                 *jac = (PC_HYPRE*)pc->data;
409:   PetscErrorCode           ierr;

412:   PCReset_HYPRE(pc);
413:   if (jac->destroy) PetscStackCall("HYPRE_DestroyXXX",(*jac->destroy)(jac->hsolver););
414:   PetscFree(jac->hypre_type);
415:   if (jac->comm_hypre != MPI_COMM_NULL) { MPI_Comm_free(&(jac->comm_hypre));}
416:   PetscFree(pc->data);

418:   PetscObjectChangeTypeName((PetscObject)pc,0);
419:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetType_C",NULL);
420:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPREGetType_C",NULL);
421:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetCoordinates_C",NULL);
422:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteGradient_C",NULL);
423:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteCurl_C",NULL);
424:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetInterpolations_C",NULL);
425:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetConstantEdgeVectors_C",NULL);
426:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetPoissonMatrix_C",NULL);
427:   return(0);
428: }

430: /* --------------------------------------------------------------------------------------------*/
431: static PetscErrorCode PCSetFromOptions_HYPRE_Pilut(PetscOptionItems *PetscOptionsObject,PC pc)
432: {
433:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
435:   PetscBool      flag;

438:   PetscOptionsHead(PetscOptionsObject,"HYPRE Pilut Options");
439:   PetscOptionsInt("-pc_hypre_pilut_maxiter","Number of iterations","None",jac->maxiter,&jac->maxiter,&flag);
440:   if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetMaxIter,(jac->hsolver,jac->maxiter));
441:   PetscOptionsReal("-pc_hypre_pilut_tol","Drop tolerance","None",jac->tol,&jac->tol,&flag);
442:   if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetDropTolerance,(jac->hsolver,jac->tol));
443:   PetscOptionsInt("-pc_hypre_pilut_factorrowsize","FactorRowSize","None",jac->factorrowsize,&jac->factorrowsize,&flag);
444:   if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetFactorRowSize,(jac->hsolver,jac->factorrowsize));
445:   PetscOptionsTail();
446:   return(0);
447: }

449: static PetscErrorCode PCView_HYPRE_Pilut(PC pc,PetscViewer viewer)
450: {
451:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
453:   PetscBool      iascii;

456:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
457:   if (iascii) {
458:     PetscViewerASCIIPrintf(viewer,"  HYPRE Pilut preconditioning\n");
459:     if (jac->maxiter != PETSC_DEFAULT) {
460:       PetscViewerASCIIPrintf(viewer,"    maximum number of iterations %d\n",jac->maxiter);
461:     } else {
462:       PetscViewerASCIIPrintf(viewer,"    default maximum number of iterations \n");
463:     }
464:     if (jac->tol != PETSC_DEFAULT) {
465:       PetscViewerASCIIPrintf(viewer,"    drop tolerance %g\n",(double)jac->tol);
466:     } else {
467:       PetscViewerASCIIPrintf(viewer,"    default drop tolerance \n");
468:     }
469:     if (jac->factorrowsize != PETSC_DEFAULT) {
470:       PetscViewerASCIIPrintf(viewer,"    factor row size %d\n",jac->factorrowsize);
471:     } else {
472:       PetscViewerASCIIPrintf(viewer,"    default factor row size \n");
473:     }
474:   }
475:   return(0);
476: }

478: /* --------------------------------------------------------------------------------------------*/

480: static PetscErrorCode PCApplyTranspose_HYPRE_BoomerAMG(PC pc,Vec b,Vec x)
481: {
482:   PC_HYPRE           *jac = (PC_HYPRE*)pc->data;
483:   Mat_HYPRE          *hjac = (Mat_HYPRE*)(jac->hpmat->data);
484:   PetscErrorCode     ierr;
485:   HYPRE_ParCSRMatrix hmat;
486:   PetscScalar        *xv;
487:   const PetscScalar  *bv;
488:   HYPRE_ParVector    jbv,jxv;
489:   PetscScalar        *sbv,*sxv;
490:   PetscInt           hierr;

493:   PetscCitationsRegister(hypreCitation,&cite);
494:   VecSet(x,0.0);
495:   VecGetArrayRead(b,&bv);
496:   VecGetArray(x,&xv);
497:   VecHYPRE_ParVectorReplacePointer(hjac->b,(PetscScalar*)bv,sbv);
498:   VecHYPRE_ParVectorReplacePointer(hjac->x,xv,sxv);

500:   PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hjac->ij,(void**)&hmat));
501:   PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->b,(void**)&jbv));
502:   PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->x,(void**)&jxv));

504:   hHYPRE_BoomerAMGSolveT(jac->hsolver,hmat,jbv,jxv);
505:   /* error code of 1 in BoomerAMG merely means convergence not achieved */
506:   if (hierr && (hierr != 1)) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"Error in HYPRE solver, error code %d",hierr);
507:   if (hierr) hypre__global_error = 0;

509:   VecHYPRE_ParVectorReplacePointer(hjac->b,sbv,bv);
510:   VecHYPRE_ParVectorReplacePointer(hjac->x,sxv,xv);
511:   VecRestoreArray(x,&xv);
512:   VecRestoreArrayRead(b,&bv);
513:   return(0);
514: }

516: /* static array length */
517: #define ALEN(a) (sizeof(a)/sizeof((a)[0]))

519: static const char *HYPREBoomerAMGCycleType[]   = {"","V","W"};
520: static const char *HYPREBoomerAMGCoarsenType[] = {"CLJP","Ruge-Stueben","","modifiedRuge-Stueben","","","Falgout", "", "PMIS", "", "HMIS"};
521: static const char *HYPREBoomerAMGMeasureType[] = {"local","global"};
522: /* The following corresponds to HYPRE_BoomerAMGSetRelaxType which has many missing numbers in the enum */
523: static const char *HYPREBoomerAMGSmoothType[]   = {"Schwarz-smoothers","Pilut","ParaSails","Euclid"};
524: static const char *HYPREBoomerAMGRelaxType[]   = {"Jacobi","sequential-Gauss-Seidel","seqboundary-Gauss-Seidel","SOR/Jacobi","backward-SOR/Jacobi",
525:                                                   "" /* [5] hybrid chaotic Gauss-Seidel (works only with OpenMP) */,"symmetric-SOR/Jacobi",
526:                                                   "" /* 7 */,"l1scaled-SOR/Jacobi","Gaussian-elimination",
527:                                                   "" /* 10 */, "" /* 11 */, "" /* 12 */, "" /* 13 */, "" /* 14 */,
528:                                                   "CG" /* non-stationary */,"Chebyshev","FCF-Jacobi","l1scaled-Jacobi"};
529: static const char *HYPREBoomerAMGInterpType[]  = {"classical", "", "", "direct", "multipass", "multipass-wts", "ext+i",
530:                                                   "ext+i-cc", "standard", "standard-wts", "block", "block-wtd", "FF", "FF1"};
531: static PetscErrorCode PCSetFromOptions_HYPRE_BoomerAMG(PetscOptionItems *PetscOptionsObject,PC pc)
532: {
533:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
535:   PetscInt       n,indx,level;
536:   PetscBool      flg, tmp_truth;
537:   double         tmpdbl, twodbl[2];

540:   PetscOptionsHead(PetscOptionsObject,"HYPRE BoomerAMG Options");
541:   PetscOptionsEList("-pc_hypre_boomeramg_cycle_type","Cycle type","None",HYPREBoomerAMGCycleType+1,2,HYPREBoomerAMGCycleType[jac->cycletype],&indx,&flg);
542:   if (flg) {
543:     jac->cycletype = indx+1;
544:     PetscStackCallStandard(HYPRE_BoomerAMGSetCycleType,(jac->hsolver,jac->cycletype));
545:   }
546:   PetscOptionsInt("-pc_hypre_boomeramg_max_levels","Number of levels (of grids) allowed","None",jac->maxlevels,&jac->maxlevels,&flg);
547:   if (flg) {
548:     if (jac->maxlevels < 2) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of levels %d must be at least two",jac->maxlevels);
549:     PetscStackCallStandard(HYPRE_BoomerAMGSetMaxLevels,(jac->hsolver,jac->maxlevels));
550:   }
551:   PetscOptionsInt("-pc_hypre_boomeramg_max_iter","Maximum iterations used PER hypre call","None",jac->maxiter,&jac->maxiter,&flg);
552:   if (flg) {
553:     if (jac->maxiter < 1) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of iterations %d must be at least one",jac->maxiter);
554:     PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,jac->maxiter));
555:   }
556:   PetscOptionsScalar("-pc_hypre_boomeramg_tol","Convergence tolerance PER hypre call (0.0 = use a fixed number of iterations)","None",jac->tol,&jac->tol,&flg);
557:   if (flg) {
558:     if (jac->tol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Tolerance %g must be greater than or equal to zero",(double)jac->tol);
559:     PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,jac->tol));
560:   }

562:   PetscOptionsScalar("-pc_hypre_boomeramg_truncfactor","Truncation factor for interpolation (0=no truncation)","None",jac->truncfactor,&jac->truncfactor,&flg);
563:   if (flg) {
564:     if (jac->truncfactor < 0.0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Truncation factor %g must be great than or equal zero",(double)jac->truncfactor);
565:     PetscStackCallStandard(HYPRE_BoomerAMGSetTruncFactor,(jac->hsolver,jac->truncfactor));
566:   }

568:   PetscOptionsInt("-pc_hypre_boomeramg_P_max","Max elements per row for interpolation operator (0=unlimited)","None",jac->pmax,&jac->pmax,&flg);
569:   if (flg) {
570:     if (jac->pmax < 0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"P_max %g must be greater than or equal to zero",(double)jac->pmax);
571:     PetscStackCallStandard(HYPRE_BoomerAMGSetPMaxElmts,(jac->hsolver,jac->pmax));
572:   }

574:   PetscOptionsInt("-pc_hypre_boomeramg_agg_nl","Number of levels of aggressive coarsening","None",jac->agg_nl,&jac->agg_nl,&flg);
575:   if (flg) {
576:     if (jac->agg_nl < 0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of levels %g must be greater than or equal to zero",(double)jac->agg_nl);

578:     PetscStackCallStandard(HYPRE_BoomerAMGSetAggNumLevels,(jac->hsolver,jac->agg_nl));
579:   }


582:   PetscOptionsInt("-pc_hypre_boomeramg_agg_num_paths","Number of paths for aggressive coarsening","None",jac->agg_num_paths,&jac->agg_num_paths,&flg);
583:   if (flg) {
584:     if (jac->agg_num_paths < 1) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of paths %g must be greater than or equal to 1",(double)jac->agg_num_paths);

586:     PetscStackCallStandard(HYPRE_BoomerAMGSetNumPaths,(jac->hsolver,jac->agg_num_paths));
587:   }


590:   PetscOptionsScalar("-pc_hypre_boomeramg_strong_threshold","Threshold for being strongly connected","None",jac->strongthreshold,&jac->strongthreshold,&flg);
591:   if (flg) {
592:     if (jac->strongthreshold < 0.0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Strong threshold %g must be great than or equal zero",(double)jac->strongthreshold);
593:     PetscStackCallStandard(HYPRE_BoomerAMGSetStrongThreshold,(jac->hsolver,jac->strongthreshold));
594:   }
595:   PetscOptionsScalar("-pc_hypre_boomeramg_max_row_sum","Maximum row sum","None",jac->maxrowsum,&jac->maxrowsum,&flg);
596:   if (flg) {
597:     if (jac->maxrowsum < 0.0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Maximum row sum %g must be greater than zero",(double)jac->maxrowsum);
598:     if (jac->maxrowsum > 1.0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Maximum row sum %g must be less than or equal one",(double)jac->maxrowsum);
599:     PetscStackCallStandard(HYPRE_BoomerAMGSetMaxRowSum,(jac->hsolver,jac->maxrowsum));
600:   }

602:   /* Grid sweeps */
603:   PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_all","Number of sweeps for the up and down grid levels","None",jac->gridsweeps[0],&indx,&flg);
604:   if (flg) {
605:     PetscStackCallStandard(HYPRE_BoomerAMGSetNumSweeps,(jac->hsolver,indx));
606:     /* modify the jac structure so we can view the updated options with PC_View */
607:     jac->gridsweeps[0] = indx;
608:     jac->gridsweeps[1] = indx;
609:     /*defaults coarse to 1 */
610:     jac->gridsweeps[2] = 1;
611:   }

613:   PetscOptionsInt("-pc_hypre_boomeramg_nodal_coarsen","Use a nodal based coarsening 1-6","HYPRE_BoomerAMGSetNodal",jac->nodal_coarsening, &jac->nodal_coarsening,&flg);
614:   if (flg) {
615:     PetscStackCallStandard(HYPRE_BoomerAMGSetNodal,(jac->hsolver,jac->nodal_coarsening));
616:   }

618:   PetscOptionsInt("-pc_hypre_boomeramg_vec_interp_variant","Variant of algorithm 1-3","HYPRE_BoomerAMGSetInterpVecVariant",jac->vec_interp_variant, &jac->vec_interp_variant,&flg);
619:   if (flg) {
620:     PetscStackCallStandard(HYPRE_BoomerAMGSetInterpVecVariant,(jac->hsolver,jac->vec_interp_variant));
621:   }

623:   PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_down","Number of sweeps for the down cycles","None",jac->gridsweeps[0], &indx,&flg);
624:   if (flg) {
625:     PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,(jac->hsolver,indx, 1));
626:     jac->gridsweeps[0] = indx;
627:   }
628:   PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_up","Number of sweeps for the up cycles","None",jac->gridsweeps[1],&indx,&flg);
629:   if (flg) {
630:     PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,(jac->hsolver,indx, 2));
631:     jac->gridsweeps[1] = indx;
632:   }
633:   PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_coarse","Number of sweeps for the coarse level","None",jac->gridsweeps[2],&indx,&flg);
634:   if (flg) {
635:     PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,(jac->hsolver,indx, 3));
636:     jac->gridsweeps[2] = indx;
637:   }

639:   /* Smooth type */
640:   PetscOptionsEList("-pc_hypre_boomeramg_smooth_type","Enable more complex smoothers","None",HYPREBoomerAMGSmoothType,ALEN(HYPREBoomerAMGSmoothType),HYPREBoomerAMGSmoothType[0],&indx,&flg);
641:   if (flg) {
642:     jac->smoothtype = indx;
643:     PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothType,(jac->hsolver,indx+6));
644:     jac->smoothnumlevels = 25;
645:     PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothNumLevels,(jac->hsolver,25));
646:   }

648:   /* Number of smoothing levels */
649:   PetscOptionsInt("-pc_hypre_boomeramg_smooth_num_levels","Number of levels on which more complex smoothers are used","None",25,&indx,&flg);
650:   if (flg && (jac->smoothtype != -1)) {
651:     jac->smoothnumlevels = indx;
652:     PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothNumLevels,(jac->hsolver,indx));
653:   }

655:   /* Number of levels for ILU(k) for Euclid */
656:   PetscOptionsInt("-pc_hypre_boomeramg_eu_level","Number of levels for ILU(k) in Euclid smoother","None",0,&indx,&flg);
657:   if (flg && (jac->smoothtype == 3)) {
658:     jac->eu_level = indx;
659:     PetscStackCallStandard(HYPRE_BoomerAMGSetEuLevel,(jac->hsolver,indx));
660:   }

662:   /* Filter for ILU(k) for Euclid */
663:   double droptolerance;
664:   PetscOptionsScalar("-pc_hypre_boomeramg_eu_droptolerance","Drop tolerance for ILU(k) in Euclid smoother","None",0,&droptolerance,&flg);
665:   if (flg && (jac->smoothtype == 3)) {
666:     jac->eu_droptolerance = droptolerance;
667:     PetscStackCallStandard(HYPRE_BoomerAMGSetEuLevel,(jac->hsolver,droptolerance));
668:   }

670:   /* Use Block Jacobi ILUT for Euclid */
671:   PetscOptionsBool("-pc_hypre_boomeramg_eu_bj", "Use Block Jacobi for ILU in Euclid smoother?", "None", PETSC_FALSE, &tmp_truth, &flg);
672:   if (flg && (jac->smoothtype == 3)) {
673:     jac->eu_bj = tmp_truth;
674:     PetscStackCallStandard(HYPRE_BoomerAMGSetEuBJ,(jac->hsolver,jac->eu_bj));
675:   }

677:   /* Relax type */
678:   PetscOptionsEList("-pc_hypre_boomeramg_relax_type_all","Relax type for the up and down cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg);
679:   if (flg) {
680:     jac->relaxtype[0] = jac->relaxtype[1]  = indx;
681:     PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxType,(jac->hsolver, indx));
682:     /* by default, coarse type set to 9 */
683:     jac->relaxtype[2] = 9;

685:   }
686:   PetscOptionsEList("-pc_hypre_boomeramg_relax_type_down","Relax type for the down cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg);
687:   if (flg) {
688:     jac->relaxtype[0] = indx;
689:     PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,(jac->hsolver, indx, 1));
690:   }
691:   PetscOptionsEList("-pc_hypre_boomeramg_relax_type_up","Relax type for the up cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg);
692:   if (flg) {
693:     jac->relaxtype[1] = indx;
694:     PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,(jac->hsolver, indx, 2));
695:   }
696:   PetscOptionsEList("-pc_hypre_boomeramg_relax_type_coarse","Relax type on coarse grid","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[9],&indx,&flg);
697:   if (flg) {
698:     jac->relaxtype[2] = indx;
699:     PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,(jac->hsolver, indx, 3));
700:   }

702:   /* Relaxation Weight */
703:   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);
704:   if (flg) {
705:     PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxWt,(jac->hsolver,tmpdbl));
706:     jac->relaxweight = tmpdbl;
707:   }

709:   n         = 2;
710:   twodbl[0] = twodbl[1] = 1.0;
711:   PetscOptionsRealArray("-pc_hypre_boomeramg_relax_weight_level","Set the relaxation weight for a particular level (weight,level)","None",twodbl, &n, &flg);
712:   if (flg) {
713:     if (n == 2) {
714:       indx =  (int)PetscAbsReal(twodbl[1]);
715:       PetscStackCallStandard(HYPRE_BoomerAMGSetLevelRelaxWt,(jac->hsolver,twodbl[0],indx));
716:     } else SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Relax weight level: you must provide 2 values separated by a comma (and no space), you provided %d",n);
717:   }

719:   /* Outer relaxation Weight */
720:   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);
721:   if (flg) {
722:     PetscStackCallStandard(HYPRE_BoomerAMGSetOuterWt,(jac->hsolver, tmpdbl));
723:     jac->outerrelaxweight = tmpdbl;
724:   }

726:   n         = 2;
727:   twodbl[0] = twodbl[1] = 1.0;
728:   PetscOptionsRealArray("-pc_hypre_boomeramg_outer_relax_weight_level","Set the outer relaxation weight for a particular level (weight,level)","None",twodbl, &n, &flg);
729:   if (flg) {
730:     if (n == 2) {
731:       indx =  (int)PetscAbsReal(twodbl[1]);
732:       PetscStackCallStandard(HYPRE_BoomerAMGSetLevelOuterWt,(jac->hsolver, twodbl[0], indx));
733:     } else SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Relax weight outer level: You must provide 2 values separated by a comma (and no space), you provided %d",n);
734:   }

736:   /* the Relax Order */
737:   PetscOptionsBool("-pc_hypre_boomeramg_no_CF", "Do not use CF-relaxation", "None", PETSC_FALSE, &tmp_truth, &flg);

739:   if (flg && tmp_truth) {
740:     jac->relaxorder = 0;
741:     PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxOrder,(jac->hsolver, jac->relaxorder));
742:   }
743:   PetscOptionsEList("-pc_hypre_boomeramg_measure_type","Measure type","None",HYPREBoomerAMGMeasureType,ALEN(HYPREBoomerAMGMeasureType),HYPREBoomerAMGMeasureType[0],&indx,&flg);
744:   if (flg) {
745:     jac->measuretype = indx;
746:     PetscStackCallStandard(HYPRE_BoomerAMGSetMeasureType,(jac->hsolver,jac->measuretype));
747:   }
748:   /* update list length 3/07 */
749:   PetscOptionsEList("-pc_hypre_boomeramg_coarsen_type","Coarsen type","None",HYPREBoomerAMGCoarsenType,ALEN(HYPREBoomerAMGCoarsenType),HYPREBoomerAMGCoarsenType[6],&indx,&flg);
750:   if (flg) {
751:     jac->coarsentype = indx;
752:     PetscStackCallStandard(HYPRE_BoomerAMGSetCoarsenType,(jac->hsolver,jac->coarsentype));
753:   }

755:   /* new 3/07 */
756:   PetscOptionsEList("-pc_hypre_boomeramg_interp_type","Interpolation type","None",HYPREBoomerAMGInterpType,ALEN(HYPREBoomerAMGInterpType),HYPREBoomerAMGInterpType[0],&indx,&flg);
757:   if (flg) {
758:     jac->interptype = indx;
759:     PetscStackCallStandard(HYPRE_BoomerAMGSetInterpType,(jac->hsolver,jac->interptype));
760:   }

762:   PetscOptionsName("-pc_hypre_boomeramg_print_statistics","Print statistics","None",&flg);
763:   if (flg) {
764:     level = 3;
765:     PetscOptionsInt("-pc_hypre_boomeramg_print_statistics","Print statistics","None",level,&level,NULL);

767:     jac->printstatistics = PETSC_TRUE;
768:     PetscStackCallStandard(HYPRE_BoomerAMGSetPrintLevel,(jac->hsolver,level));
769:   }

771:   PetscOptionsName("-pc_hypre_boomeramg_print_debug","Print debug information","None",&flg);
772:   if (flg) {
773:     level = 3;
774:     PetscOptionsInt("-pc_hypre_boomeramg_print_debug","Print debug information","None",level,&level,NULL);

776:     jac->printstatistics = PETSC_TRUE;
777:     PetscStackCallStandard(HYPRE_BoomerAMGSetDebugFlag,(jac->hsolver,level));
778:   }

780:   PetscOptionsBool("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None", PETSC_FALSE, &tmp_truth, &flg);
781:   if (flg && tmp_truth) {
782:     PetscInt tmp_int;
783:     PetscOptionsInt("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None",jac->nodal_relax_levels,&tmp_int,&flg);
784:     if (flg) jac->nodal_relax_levels = tmp_int;
785:     PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothType,(jac->hsolver,6));
786:     PetscStackCallStandard(HYPRE_BoomerAMGSetDomainType,(jac->hsolver,1));
787:     PetscStackCallStandard(HYPRE_BoomerAMGSetOverlap,(jac->hsolver,0));
788:     PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothNumLevels,(jac->hsolver,jac->nodal_relax_levels));
789:   }

791:   PetscOptionsTail();
792:   return(0);
793: }

795: static PetscErrorCode PCApplyRichardson_HYPRE_BoomerAMG(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its,PetscBool guesszero,PetscInt *outits,PCRichardsonConvergedReason *reason)
796: {
797:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
799:   PetscInt       oits;

802:   PetscCitationsRegister(hypreCitation,&cite);
803:   PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,its*jac->maxiter));
804:   PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,rtol));
805:   jac->applyrichardson = PETSC_TRUE;
806:   PCApply_HYPRE(pc,b,y);
807:   jac->applyrichardson = PETSC_FALSE;
808:   PetscStackCallStandard(HYPRE_BoomerAMGGetNumIterations,(jac->hsolver,(HYPRE_Int *)&oits));
809:   *outits = oits;
810:   if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
811:   else             *reason = PCRICHARDSON_CONVERGED_RTOL;
812:   PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,jac->tol));
813:   PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,jac->maxiter));
814:   return(0);
815: }


818: static PetscErrorCode PCView_HYPRE_BoomerAMG(PC pc,PetscViewer viewer)
819: {
820:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
822:   PetscBool      iascii;

825:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
826:   if (iascii) {
827:     PetscViewerASCIIPrintf(viewer,"  HYPRE BoomerAMG preconditioning\n");
828:     PetscViewerASCIIPrintf(viewer,"    Cycle type %s\n",HYPREBoomerAMGCycleType[jac->cycletype]);
829:     PetscViewerASCIIPrintf(viewer,"    Maximum number of levels %d\n",jac->maxlevels);
830:     PetscViewerASCIIPrintf(viewer,"    Maximum number of iterations PER hypre call %d\n",jac->maxiter);
831:     PetscViewerASCIIPrintf(viewer,"    Convergence tolerance PER hypre call %g\n",(double)jac->tol);
832:     PetscViewerASCIIPrintf(viewer,"    Threshold for strong coupling %g\n",(double)jac->strongthreshold);
833:     PetscViewerASCIIPrintf(viewer,"    Interpolation truncation factor %g\n",(double)jac->truncfactor);
834:     PetscViewerASCIIPrintf(viewer,"    Interpolation: max elements per row %d\n",jac->pmax);
835:     PetscViewerASCIIPrintf(viewer,"    Number of levels of aggressive coarsening %d\n",jac->agg_nl);
836:     PetscViewerASCIIPrintf(viewer,"    Number of paths for aggressive coarsening %d\n",jac->agg_num_paths);

838:     PetscViewerASCIIPrintf(viewer,"    Maximum row sums %g\n",(double)jac->maxrowsum);

840:     PetscViewerASCIIPrintf(viewer,"    Sweeps down         %d\n",jac->gridsweeps[0]);
841:     PetscViewerASCIIPrintf(viewer,"    Sweeps up           %d\n",jac->gridsweeps[1]);
842:     PetscViewerASCIIPrintf(viewer,"    Sweeps on coarse    %d\n",jac->gridsweeps[2]);

844:     PetscViewerASCIIPrintf(viewer,"    Relax down          %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[0]]);
845:     PetscViewerASCIIPrintf(viewer,"    Relax up            %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[1]]);
846:     PetscViewerASCIIPrintf(viewer,"    Relax on coarse     %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[2]]);

848:     PetscViewerASCIIPrintf(viewer,"    Relax weight  (all)      %g\n",(double)jac->relaxweight);
849:     PetscViewerASCIIPrintf(viewer,"    Outer relax weight (all) %g\n",(double)jac->outerrelaxweight);

851:     if (jac->relaxorder) {
852:       PetscViewerASCIIPrintf(viewer,"    Using CF-relaxation\n");
853:     } else {
854:       PetscViewerASCIIPrintf(viewer,"    Not using CF-relaxation\n");
855:     }
856:     if (jac->smoothtype!=-1) {
857:       PetscViewerASCIIPrintf(viewer,"    Smooth type          %s\n",HYPREBoomerAMGSmoothType[jac->smoothtype]);
858:       PetscViewerASCIIPrintf(viewer,"    Smooth num levels    %d\n",jac->smoothnumlevels);
859:     } else {
860:       PetscViewerASCIIPrintf(viewer,"    Not using more complex smoothers.\n");
861:     }
862:     if (jac->smoothtype==3) {
863:       PetscViewerASCIIPrintf(viewer,"    Euclid ILU(k) levels %d\n",jac->eu_level);
864:       PetscViewerASCIIPrintf(viewer,"    Euclid ILU(k) drop tolerance %g\n",jac->eu_droptolerance);
865:       PetscViewerASCIIPrintf(viewer,"    Euclid ILU use Block-Jacobi? %d\n",jac->eu_bj);
866:     }
867:     PetscViewerASCIIPrintf(viewer,"    Measure type        %s\n",HYPREBoomerAMGMeasureType[jac->measuretype]);
868:     PetscViewerASCIIPrintf(viewer,"    Coarsen type        %s\n",HYPREBoomerAMGCoarsenType[jac->coarsentype]);
869:     PetscViewerASCIIPrintf(viewer,"    Interpolation type  %s\n",HYPREBoomerAMGInterpType[jac->interptype]);
870:     if (jac->nodal_coarsening) {
871:       PetscViewerASCIIPrintf(viewer,"    Using nodal coarsening (with HYPRE_BOOMERAMGSetNodal() %D\n",jac->nodal_coarsening);
872:     }
873:     if (jac->vec_interp_variant) {
874:       PetscViewerASCIIPrintf(viewer,"    HYPRE_BoomerAMGSetInterpVecVariant() %D\n",jac->vec_interp_variant);
875:     }
876:     if (jac->nodal_relax) {
877:       PetscViewerASCIIPrintf(viewer,"    Using nodal relaxation via Schwarz smoothing on levels %d\n",jac->nodal_relax_levels);
878:     }
879:   }
880:   return(0);
881: }

883: /* --------------------------------------------------------------------------------------------*/
884: static PetscErrorCode PCSetFromOptions_HYPRE_ParaSails(PetscOptionItems *PetscOptionsObject,PC pc)
885: {
886:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
888:   PetscInt       indx;
889:   PetscBool      flag;
890:   const char     *symtlist[] = {"nonsymmetric","SPD","nonsymmetric,SPD"};

893:   PetscOptionsHead(PetscOptionsObject,"HYPRE ParaSails Options");
894:   PetscOptionsInt("-pc_hypre_parasails_nlevels","Number of number of levels","None",jac->nlevels,&jac->nlevels,0);
895:   PetscOptionsReal("-pc_hypre_parasails_thresh","Threshold","None",jac->threshhold,&jac->threshhold,&flag);
896:   if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetParams,(jac->hsolver,jac->threshhold,jac->nlevels));

898:   PetscOptionsReal("-pc_hypre_parasails_filter","filter","None",jac->filter,&jac->filter,&flag);
899:   if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetFilter,(jac->hsolver,jac->filter));

901:   PetscOptionsReal("-pc_hypre_parasails_loadbal","Load balance","None",jac->loadbal,&jac->loadbal,&flag);
902:   if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetLoadbal,(jac->hsolver,jac->loadbal));

904:   PetscOptionsBool("-pc_hypre_parasails_logging","Print info to screen","None",(PetscBool)jac->logging,(PetscBool*)&jac->logging,&flag);
905:   if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetLogging,(jac->hsolver,jac->logging));

907:   PetscOptionsBool("-pc_hypre_parasails_reuse","Reuse nonzero pattern in preconditioner","None",(PetscBool)jac->ruse,(PetscBool*)&jac->ruse,&flag);
908:   if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetReuse,(jac->hsolver,jac->ruse));

910:   PetscOptionsEList("-pc_hypre_parasails_sym","Symmetry of matrix and preconditioner","None",symtlist,ALEN(symtlist),symtlist[0],&indx,&flag);
911:   if (flag) {
912:     jac->symt = indx;
913:     PetscStackCallStandard(HYPRE_ParaSailsSetSym,(jac->hsolver,jac->symt));
914:   }

916:   PetscOptionsTail();
917:   return(0);
918: }

920: static PetscErrorCode PCView_HYPRE_ParaSails(PC pc,PetscViewer viewer)
921: {
922:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
924:   PetscBool      iascii;
925:   const char     *symt = 0;;

928:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
929:   if (iascii) {
930:     PetscViewerASCIIPrintf(viewer,"  HYPRE ParaSails preconditioning\n");
931:     PetscViewerASCIIPrintf(viewer,"    nlevels %d\n",jac->nlevels);
932:     PetscViewerASCIIPrintf(viewer,"    threshold %g\n",(double)jac->threshhold);
933:     PetscViewerASCIIPrintf(viewer,"    filter %g\n",(double)jac->filter);
934:     PetscViewerASCIIPrintf(viewer,"    load balance %g\n",(double)jac->loadbal);
935:     PetscViewerASCIIPrintf(viewer,"    reuse nonzero structure %s\n",PetscBools[jac->ruse]);
936:     PetscViewerASCIIPrintf(viewer,"    print info to screen %s\n",PetscBools[jac->logging]);
937:     if (!jac->symt) symt = "nonsymmetric matrix and preconditioner";
938:     else if (jac->symt == 1) symt = "SPD matrix and preconditioner";
939:     else if (jac->symt == 2) symt = "nonsymmetric matrix but SPD preconditioner";
940:     else SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG,"Unknown HYPRE ParaSails symmetric option %d",jac->symt);
941:     PetscViewerASCIIPrintf(viewer,"    %s\n",symt);
942:   }
943:   return(0);
944: }
945: /* --------------------------------------------------------------------------------------------*/
946: static PetscErrorCode PCSetFromOptions_HYPRE_AMS(PetscOptionItems *PetscOptionsObject,PC pc)
947: {
948:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
950:   PetscInt       n;
951:   PetscBool      flag,flag2,flag3,flag4;

954:   PetscOptionsHead(PetscOptionsObject,"HYPRE AMS Options");
955:   PetscOptionsInt("-pc_hypre_ams_print_level","Debugging output level for AMS","None",jac->as_print,&jac->as_print,&flag);
956:   if (flag) PetscStackCallStandard(HYPRE_AMSSetPrintLevel,(jac->hsolver,jac->as_print));
957:   PetscOptionsInt("-pc_hypre_ams_max_iter","Maximum number of AMS multigrid iterations within PCApply","None",jac->as_max_iter,&jac->as_max_iter,&flag);
958:   if (flag) PetscStackCallStandard(HYPRE_AMSSetMaxIter,(jac->hsolver,jac->as_max_iter));
959:   PetscOptionsInt("-pc_hypre_ams_cycle_type","Cycle type for AMS multigrid","None",jac->ams_cycle_type,&jac->ams_cycle_type,&flag);
960:   if (flag) PetscStackCallStandard(HYPRE_AMSSetCycleType,(jac->hsolver,jac->ams_cycle_type));
961:   PetscOptionsReal("-pc_hypre_ams_tol","Error tolerance for AMS multigrid","None",jac->as_tol,&jac->as_tol,&flag);
962:   if (flag) PetscStackCallStandard(HYPRE_AMSSetTol,(jac->hsolver,jac->as_tol));
963:   PetscOptionsInt("-pc_hypre_ams_relax_type","Relaxation type for AMS smoother","None",jac->as_relax_type,&jac->as_relax_type,&flag);
964:   PetscOptionsInt("-pc_hypre_ams_relax_times","Number of relaxation steps for AMS smoother","None",jac->as_relax_times,&jac->as_relax_times,&flag2);
965:   PetscOptionsReal("-pc_hypre_ams_relax_weight","Relaxation weight for AMS smoother","None",jac->as_relax_weight,&jac->as_relax_weight,&flag3);
966:   PetscOptionsReal("-pc_hypre_ams_omega","SSOR coefficient for AMS smoother","None",jac->as_omega,&jac->as_omega,&flag4);
967:   if (flag || flag2 || flag3 || flag4) {
968:     PetscStackCallStandard(HYPRE_AMSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type,
969:                                                                       jac->as_relax_times,
970:                                                                       jac->as_relax_weight,
971:                                                                       jac->as_omega));
972:   }
973:   PetscOptionsReal("-pc_hypre_ams_amg_alpha_theta","Threshold for strong coupling of vector Poisson AMG solver","None",jac->as_amg_alpha_theta,&jac->as_amg_alpha_theta,&flag);
974:   n = 5;
975:   PetscOptionsIntArray("-pc_hypre_ams_amg_alpha_options","AMG options for vector Poisson","None",jac->as_amg_alpha_opts,&n,&flag2);
976:   if (flag || flag2) {
977:     PetscStackCallStandard(HYPRE_AMSSetAlphaAMGOptions,(jac->hsolver,jac->as_amg_alpha_opts[0],       /* AMG coarsen type */
978:                                                                      jac->as_amg_alpha_opts[1],       /* AMG agg_levels */
979:                                                                      jac->as_amg_alpha_opts[2],       /* AMG relax_type */
980:                                                                      jac->as_amg_alpha_theta,
981:                                                                      jac->as_amg_alpha_opts[3],       /* AMG interp_type */
982:                                                                      jac->as_amg_alpha_opts[4]));     /* AMG Pmax */
983:   }
984:   PetscOptionsReal("-pc_hypre_ams_amg_beta_theta","Threshold for strong coupling of scalar Poisson AMG solver","None",jac->as_amg_beta_theta,&jac->as_amg_beta_theta,&flag);
985:   n = 5;
986:   PetscOptionsIntArray("-pc_hypre_ams_amg_beta_options","AMG options for scalar Poisson solver","None",jac->as_amg_beta_opts,&n,&flag2);
987:   if (flag || flag2) {
988:     PetscStackCallStandard(HYPRE_AMSSetBetaAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0],       /* AMG coarsen type */
989:                                                                     jac->as_amg_beta_opts[1],       /* AMG agg_levels */
990:                                                                     jac->as_amg_beta_opts[2],       /* AMG relax_type */
991:                                                                     jac->as_amg_beta_theta,
992:                                                                     jac->as_amg_beta_opts[3],       /* AMG interp_type */
993:                                                                     jac->as_amg_beta_opts[4]));     /* AMG Pmax */
994:   }
995:   PetscOptionsInt("-pc_hypre_ams_projection_frequency","Frequency at which a projection onto the compatible subspace for problems with zero conductivity regions is performed","None",jac->ams_proj_freq,&jac->ams_proj_freq,&flag);
996:   if (flag) { /* override HYPRE's default only if the options is used */
997:     PetscStackCallStandard(HYPRE_AMSSetProjectionFrequency,(jac->hsolver,jac->ams_proj_freq));
998:   }
999:   PetscOptionsTail();
1000:   return(0);
1001: }

1003: static PetscErrorCode PCView_HYPRE_AMS(PC pc,PetscViewer viewer)
1004: {
1005:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
1007:   PetscBool      iascii;

1010:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
1011:   if (iascii) {
1012:     PetscViewerASCIIPrintf(viewer,"  HYPRE AMS preconditioning\n");
1013:     PetscViewerASCIIPrintf(viewer,"    subspace iterations per application %d\n",jac->as_max_iter);
1014:     PetscViewerASCIIPrintf(viewer,"    subspace cycle type %d\n",jac->ams_cycle_type);
1015:     PetscViewerASCIIPrintf(viewer,"    subspace iteration tolerance %g\n",jac->as_tol);
1016:     PetscViewerASCIIPrintf(viewer,"    smoother type %d\n",jac->as_relax_type);
1017:     PetscViewerASCIIPrintf(viewer,"    number of smoothing steps %d\n",jac->as_relax_times);
1018:     PetscViewerASCIIPrintf(viewer,"    smoother weight %g\n",jac->as_relax_weight);
1019:     PetscViewerASCIIPrintf(viewer,"    smoother omega %g\n",jac->as_omega);
1020:     if (jac->alpha_Poisson) {
1021:       PetscViewerASCIIPrintf(viewer,"    vector Poisson solver (passed in by user)\n");
1022:     } else {
1023:       PetscViewerASCIIPrintf(viewer,"    vector Poisson solver (computed) \n");
1024:     }
1025:     PetscViewerASCIIPrintf(viewer,"        boomerAMG coarsening type %d\n",jac->as_amg_alpha_opts[0]);
1026:     PetscViewerASCIIPrintf(viewer,"        boomerAMG levels of aggressive coarsening %d\n",jac->as_amg_alpha_opts[1]);
1027:     PetscViewerASCIIPrintf(viewer,"        boomerAMG relaxation type %d\n",jac->as_amg_alpha_opts[2]);
1028:     PetscViewerASCIIPrintf(viewer,"        boomerAMG interpolation type %d\n",jac->as_amg_alpha_opts[3]);
1029:     PetscViewerASCIIPrintf(viewer,"        boomerAMG max nonzero elements in interpolation rows %d\n",jac->as_amg_alpha_opts[4]);
1030:     PetscViewerASCIIPrintf(viewer,"        boomerAMG strength threshold %g\n",jac->as_amg_alpha_theta);
1031:     if (!jac->ams_beta_is_zero) {
1032:       if (jac->beta_Poisson) {
1033:         PetscViewerASCIIPrintf(viewer,"    scalar Poisson solver (passed in by user)\n");
1034:       } else {
1035:         PetscViewerASCIIPrintf(viewer,"    scalar Poisson solver (computed) \n");
1036:       }
1037:       PetscViewerASCIIPrintf(viewer,"        boomerAMG coarsening type %d\n",jac->as_amg_beta_opts[0]);
1038:       PetscViewerASCIIPrintf(viewer,"        boomerAMG levels of aggressive coarsening %d\n",jac->as_amg_beta_opts[1]);
1039:       PetscViewerASCIIPrintf(viewer,"        boomerAMG relaxation type %d\n",jac->as_amg_beta_opts[2]);
1040:       PetscViewerASCIIPrintf(viewer,"        boomerAMG interpolation type %d\n",jac->as_amg_beta_opts[3]);
1041:       PetscViewerASCIIPrintf(viewer,"        boomerAMG max nonzero elements in interpolation rows %d\n",jac->as_amg_beta_opts[4]);
1042:       PetscViewerASCIIPrintf(viewer,"        boomerAMG strength threshold %g\n",jac->as_amg_beta_theta);
1043:       if (jac->ams_beta_is_zero_part) {
1044:         PetscViewerASCIIPrintf(viewer,"        compatible subspace projection frequency %d (-1 HYPRE uses default)\n",jac->ams_proj_freq);
1045:       }
1046:     } else {
1047:       PetscViewerASCIIPrintf(viewer,"    scalar Poisson solver not used (zero-conductivity everywhere) \n");
1048:     }
1049:   }
1050:   return(0);
1051: }

1053: static PetscErrorCode PCSetFromOptions_HYPRE_ADS(PetscOptionItems *PetscOptionsObject,PC pc)
1054: {
1055:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
1057:   PetscInt       n;
1058:   PetscBool      flag,flag2,flag3,flag4;

1061:   PetscOptionsHead(PetscOptionsObject,"HYPRE ADS Options");
1062:   PetscOptionsInt("-pc_hypre_ads_print_level","Debugging output level for ADS","None",jac->as_print,&jac->as_print,&flag);
1063:   if (flag) PetscStackCallStandard(HYPRE_ADSSetPrintLevel,(jac->hsolver,jac->as_print));
1064:   PetscOptionsInt("-pc_hypre_ads_max_iter","Maximum number of ADS multigrid iterations within PCApply","None",jac->as_max_iter,&jac->as_max_iter,&flag);
1065:   if (flag) PetscStackCallStandard(HYPRE_ADSSetMaxIter,(jac->hsolver,jac->as_max_iter));
1066:   PetscOptionsInt("-pc_hypre_ads_cycle_type","Cycle type for ADS multigrid","None",jac->ads_cycle_type,&jac->ads_cycle_type,&flag);
1067:   if (flag) PetscStackCallStandard(HYPRE_ADSSetCycleType,(jac->hsolver,jac->ads_cycle_type));
1068:   PetscOptionsReal("-pc_hypre_ads_tol","Error tolerance for ADS multigrid","None",jac->as_tol,&jac->as_tol,&flag);
1069:   if (flag) PetscStackCallStandard(HYPRE_ADSSetTol,(jac->hsolver,jac->as_tol));
1070:   PetscOptionsInt("-pc_hypre_ads_relax_type","Relaxation type for ADS smoother","None",jac->as_relax_type,&jac->as_relax_type,&flag);
1071:   PetscOptionsInt("-pc_hypre_ads_relax_times","Number of relaxation steps for ADS smoother","None",jac->as_relax_times,&jac->as_relax_times,&flag2);
1072:   PetscOptionsReal("-pc_hypre_ads_relax_weight","Relaxation weight for ADS smoother","None",jac->as_relax_weight,&jac->as_relax_weight,&flag3);
1073:   PetscOptionsReal("-pc_hypre_ads_omega","SSOR coefficient for ADS smoother","None",jac->as_omega,&jac->as_omega,&flag4);
1074:   if (flag || flag2 || flag3 || flag4) {
1075:     PetscStackCallStandard(HYPRE_ADSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type,
1076:                                                                       jac->as_relax_times,
1077:                                                                       jac->as_relax_weight,
1078:                                                                       jac->as_omega));
1079:   }
1080:   PetscOptionsReal("-pc_hypre_ads_ams_theta","Threshold for strong coupling of AMS solver inside ADS","None",jac->as_amg_alpha_theta,&jac->as_amg_alpha_theta,&flag);
1081:   n = 5;
1082:   PetscOptionsIntArray("-pc_hypre_ads_ams_options","AMG options for AMS solver inside ADS","None",jac->as_amg_alpha_opts,&n,&flag2);
1083:   PetscOptionsInt("-pc_hypre_ads_ams_cycle_type","Cycle type for AMS solver inside ADS","None",jac->ams_cycle_type,&jac->ams_cycle_type,&flag3);
1084:   if (flag || flag2 || flag3) {
1085:     PetscStackCallStandard(HYPRE_ADSSetAMSOptions,(jac->hsolver,jac->ams_cycle_type,             /* AMS cycle type */
1086:                                                                 jac->as_amg_alpha_opts[0],       /* AMG coarsen type */
1087:                                                                 jac->as_amg_alpha_opts[1],       /* AMG agg_levels */
1088:                                                                 jac->as_amg_alpha_opts[2],       /* AMG relax_type */
1089:                                                                 jac->as_amg_alpha_theta,
1090:                                                                 jac->as_amg_alpha_opts[3],       /* AMG interp_type */
1091:                                                                 jac->as_amg_alpha_opts[4]));     /* AMG Pmax */
1092:   }
1093:   PetscOptionsReal("-pc_hypre_ads_amg_theta","Threshold for strong coupling of vector AMG solver inside ADS","None",jac->as_amg_beta_theta,&jac->as_amg_beta_theta,&flag);
1094:   n = 5;
1095:   PetscOptionsIntArray("-pc_hypre_ads_amg_options","AMG options for vector AMG solver inside ADS","None",jac->as_amg_beta_opts,&n,&flag2);
1096:   if (flag || flag2) {
1097:     PetscStackCallStandard(HYPRE_ADSSetAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0],       /* AMG coarsen type */
1098:                                                                 jac->as_amg_beta_opts[1],       /* AMG agg_levels */
1099:                                                                 jac->as_amg_beta_opts[2],       /* AMG relax_type */
1100:                                                                 jac->as_amg_beta_theta,
1101:                                                                 jac->as_amg_beta_opts[3],       /* AMG interp_type */
1102:                                                                 jac->as_amg_beta_opts[4]));     /* AMG Pmax */
1103:   }
1104:   PetscOptionsTail();
1105:   return(0);
1106: }

1108: static PetscErrorCode PCView_HYPRE_ADS(PC pc,PetscViewer viewer)
1109: {
1110:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
1112:   PetscBool      iascii;

1115:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
1116:   if (iascii) {
1117:     PetscViewerASCIIPrintf(viewer,"  HYPRE ADS preconditioning\n");
1118:     PetscViewerASCIIPrintf(viewer,"    subspace iterations per application %d\n",jac->as_max_iter);
1119:     PetscViewerASCIIPrintf(viewer,"    subspace cycle type %d\n",jac->ads_cycle_type);
1120:     PetscViewerASCIIPrintf(viewer,"    subspace iteration tolerance %g\n",jac->as_tol);
1121:     PetscViewerASCIIPrintf(viewer,"    smoother type %d\n",jac->as_relax_type);
1122:     PetscViewerASCIIPrintf(viewer,"    number of smoothing steps %d\n",jac->as_relax_times);
1123:     PetscViewerASCIIPrintf(viewer,"    smoother weight %g\n",jac->as_relax_weight);
1124:     PetscViewerASCIIPrintf(viewer,"    smoother omega %g\n",jac->as_omega);
1125:     PetscViewerASCIIPrintf(viewer,"    AMS solver using boomerAMG\n");
1126:     PetscViewerASCIIPrintf(viewer,"        subspace cycle type %d\n",jac->ams_cycle_type);
1127:     PetscViewerASCIIPrintf(viewer,"        coarsening type %d\n",jac->as_amg_alpha_opts[0]);
1128:     PetscViewerASCIIPrintf(viewer,"        levels of aggressive coarsening %d\n",jac->as_amg_alpha_opts[1]);
1129:     PetscViewerASCIIPrintf(viewer,"        relaxation type %d\n",jac->as_amg_alpha_opts[2]);
1130:     PetscViewerASCIIPrintf(viewer,"        interpolation type %d\n",jac->as_amg_alpha_opts[3]);
1131:     PetscViewerASCIIPrintf(viewer,"        max nonzero elements in interpolation rows %d\n",jac->as_amg_alpha_opts[4]);
1132:     PetscViewerASCIIPrintf(viewer,"        strength threshold %g\n",jac->as_amg_alpha_theta);
1133:     PetscViewerASCIIPrintf(viewer,"    vector Poisson solver using boomerAMG\n");
1134:     PetscViewerASCIIPrintf(viewer,"        coarsening type %d\n",jac->as_amg_beta_opts[0]);
1135:     PetscViewerASCIIPrintf(viewer,"        levels of aggressive coarsening %d\n",jac->as_amg_beta_opts[1]);
1136:     PetscViewerASCIIPrintf(viewer,"        relaxation type %d\n",jac->as_amg_beta_opts[2]);
1137:     PetscViewerASCIIPrintf(viewer,"        interpolation type %d\n",jac->as_amg_beta_opts[3]);
1138:     PetscViewerASCIIPrintf(viewer,"        max nonzero elements in interpolation rows %d\n",jac->as_amg_beta_opts[4]);
1139:     PetscViewerASCIIPrintf(viewer,"        strength threshold %g\n",jac->as_amg_beta_theta);
1140:   }
1141:   return(0);
1142: }

1144: static PetscErrorCode PCHYPRESetDiscreteGradient_HYPRE(PC pc, Mat G)
1145: {
1146:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
1147:   PetscBool      ishypre;

1151:   PetscObjectTypeCompare((PetscObject)G,MATHYPRE,&ishypre);
1152:   if (ishypre) {
1153:     PetscObjectReference((PetscObject)G);
1154:     MatDestroy(&jac->G);
1155:     jac->G = G;
1156:   } else {
1157:     MatDestroy(&jac->G);
1158:     MatConvert(G,MATHYPRE,MAT_INITIAL_MATRIX,&jac->G);
1159:   }
1160:   return(0);
1161: }

1163: /*@
1164:  PCHYPRESetDiscreteGradient - Set discrete gradient matrix

1166:    Collective on PC

1168:    Input Parameters:
1169: +  pc - the preconditioning context
1170: -  G - the discrete gradient

1172:    Level: intermediate

1174:    Notes: G should have as many rows as the number of edges and as many columns as the number of vertices in the mesh
1175:           Each row of G has 2 nonzeros, with column indexes being the global indexes of edge's endpoints: matrix entries are +1 and -1 depending on edge orientation

1177: .seealso:
1178: @*/
1179: PetscErrorCode PCHYPRESetDiscreteGradient(PC pc, Mat G)
1180: {

1187:   PetscTryMethod(pc,"PCHYPRESetDiscreteGradient_C",(PC,Mat),(pc,G));
1188:   return(0);
1189: }

1191: static PetscErrorCode PCHYPRESetDiscreteCurl_HYPRE(PC pc, Mat C)
1192: {
1193:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
1194:   PetscBool      ishypre;

1198:   PetscObjectTypeCompare((PetscObject)C,MATHYPRE,&ishypre);
1199:   if (ishypre) {
1200:     PetscObjectReference((PetscObject)C);
1201:     MatDestroy(&jac->C);
1202:     jac->C = C;
1203:   } else {
1204:     MatDestroy(&jac->C);
1205:     MatConvert(C,MATHYPRE,MAT_INITIAL_MATRIX,&jac->C);
1206:   }
1207:   return(0);
1208: }

1210: /*@
1211:  PCHYPRESetDiscreteCurl - Set discrete curl matrix

1213:    Collective on PC

1215:    Input Parameters:
1216: +  pc - the preconditioning context
1217: -  C - the discrete curl

1219:    Level: intermediate

1221:    Notes: C should have as many rows as the number of faces and as many columns as the number of edges in the mesh
1222:           Each row of G has as many nonzeros as the number of edges of a face, with column indexes being the global indexes of the corresponding edge: matrix entries are +1 and -1 depending on edge orientation with respect to the face orientation

1224: .seealso:
1225: @*/
1226: PetscErrorCode PCHYPRESetDiscreteCurl(PC pc, Mat C)
1227: {

1234:   PetscTryMethod(pc,"PCHYPRESetDiscreteCurl_C",(PC,Mat),(pc,C));
1235:   return(0);
1236: }

1238: static PetscErrorCode PCHYPRESetInterpolations_HYPRE(PC pc, PetscInt dim, Mat RT_PiFull, Mat RT_Pi[], Mat ND_PiFull, Mat ND_Pi[])
1239: {
1240:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
1241:   PetscBool      ishypre;
1243:   PetscInt       i;

1246:   MatDestroy(&jac->RT_PiFull);
1247:   MatDestroy(&jac->ND_PiFull);
1248:   for (i=0;i<3;++i) {
1249:     MatDestroy(&jac->RT_Pi[i]);
1250:     MatDestroy(&jac->ND_Pi[i]);
1251:   }

1253:   jac->dim = dim;
1254:   if (RT_PiFull) {
1255:     PetscObjectTypeCompare((PetscObject)RT_PiFull,MATHYPRE,&ishypre);
1256:     if (ishypre) {
1257:       PetscObjectReference((PetscObject)RT_PiFull);
1258:       jac->RT_PiFull = RT_PiFull;
1259:     } else {
1260:       MatConvert(RT_PiFull,MATHYPRE,MAT_INITIAL_MATRIX,&jac->RT_PiFull);
1261:     }
1262:   }
1263:   if (RT_Pi) {
1264:     for (i=0;i<dim;++i) {
1265:       if (RT_Pi[i]) {
1266:         PetscObjectTypeCompare((PetscObject)RT_Pi[i],MATHYPRE,&ishypre);
1267:         if (ishypre) {
1268:           PetscObjectReference((PetscObject)RT_Pi[i]);
1269:           jac->RT_Pi[i] = RT_Pi[i];
1270:         } else {
1271:           MatConvert(RT_Pi[i],MATHYPRE,MAT_INITIAL_MATRIX,&jac->RT_Pi[i]);
1272:         }
1273:       }
1274:     }
1275:   }
1276:   if (ND_PiFull) {
1277:     PetscObjectTypeCompare((PetscObject)ND_PiFull,MATHYPRE,&ishypre);
1278:     if (ishypre) {
1279:       PetscObjectReference((PetscObject)ND_PiFull);
1280:       jac->ND_PiFull = ND_PiFull;
1281:     } else {
1282:       MatConvert(ND_PiFull,MATHYPRE,MAT_INITIAL_MATRIX,&jac->ND_PiFull);
1283:     }
1284:   }
1285:   if (ND_Pi) {
1286:     for (i=0;i<dim;++i) {
1287:       if (ND_Pi[i]) {
1288:         PetscObjectTypeCompare((PetscObject)ND_Pi[i],MATHYPRE,&ishypre);
1289:         if (ishypre) {
1290:           PetscObjectReference((PetscObject)ND_Pi[i]);
1291:           jac->ND_Pi[i] = ND_Pi[i];
1292:         } else {
1293:           MatConvert(ND_Pi[i],MATHYPRE,MAT_INITIAL_MATRIX,&jac->ND_Pi[i]);
1294:         }
1295:       }
1296:     }
1297:   }

1299:   return(0);
1300: }

1302: /*@
1303:  PCHYPRESetInterpolations - Set interpolation matrices for AMS/ADS preconditioner

1305:    Collective on PC

1307:    Input Parameters:
1308: +  pc - the preconditioning context
1309: -  dim - the dimension of the problem, only used in AMS
1310: -  RT_PiFull - Raviart-Thomas interpolation matrix
1311: -  RT_Pi - x/y/z component of Raviart-Thomas interpolation matrix
1312: -  ND_PiFull - Nedelec interpolation matrix
1313: -  ND_Pi - x/y/z component of Nedelec interpolation matrix

1315:    Notes: For AMS, only Nedelec interpolation matrices are needed, the Raviart-Thomas interpolation matrices can be set to NULL.
1316:           For ADS, both type of interpolation matrices are needed.
1317:    Level: intermediate

1319: .seealso:
1320: @*/
1321: PetscErrorCode PCHYPRESetInterpolations(PC pc, PetscInt dim, Mat RT_PiFull, Mat RT_Pi[], Mat ND_PiFull, Mat ND_Pi[])
1322: {
1324:   PetscInt       i;

1328:   if (RT_PiFull) {
1331:   }
1332:   if (RT_Pi) {
1334:     for (i=0;i<dim;++i) {
1335:       if (RT_Pi[i]) {
1338:       }
1339:     }
1340:   }
1341:   if (ND_PiFull) {
1344:   }
1345:   if (ND_Pi) {
1347:     for (i=0;i<dim;++i) {
1348:       if (ND_Pi[i]) {
1351:       }
1352:     }
1353:   }
1354:   PetscTryMethod(pc,"PCHYPRESetInterpolations_C",(PC,PetscInt,Mat,Mat[],Mat,Mat[]),(pc,dim,RT_PiFull,RT_Pi,ND_PiFull,ND_Pi));
1355:   return(0);
1356: }

1358: static PetscErrorCode PCHYPRESetPoissonMatrix_HYPRE(PC pc, Mat A, PetscBool isalpha)
1359: {
1360:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
1361:   PetscBool      ishypre;

1365:   PetscObjectTypeCompare((PetscObject)A,MATHYPRE,&ishypre);
1366:   if (ishypre) {
1367:     if (isalpha) {
1368:       PetscObjectReference((PetscObject)A);
1369:       MatDestroy(&jac->alpha_Poisson);
1370:       jac->alpha_Poisson = A;
1371:     } else {
1372:       if (A) {
1373:         PetscObjectReference((PetscObject)A);
1374:       } else {
1375:         jac->ams_beta_is_zero = PETSC_TRUE;
1376:       }
1377:       MatDestroy(&jac->beta_Poisson);
1378:       jac->beta_Poisson = A;
1379:     }
1380:   } else {
1381:     if (isalpha) {
1382:       MatDestroy(&jac->alpha_Poisson);
1383:       MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&jac->alpha_Poisson);
1384:     } else {
1385:       if (A) {
1386:         MatDestroy(&jac->beta_Poisson);
1387:         MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&jac->beta_Poisson);
1388:       } else {
1389:         MatDestroy(&jac->beta_Poisson);
1390:         jac->ams_beta_is_zero = PETSC_TRUE;
1391:       }
1392:     }
1393:   }
1394:   return(0);
1395: }

1397: /*@
1398:  PCHYPRESetAlphaPoissonMatrix - Set vector Poisson matrix

1400:    Collective on PC

1402:    Input Parameters:
1403: +  pc - the preconditioning context
1404: -  A - the matrix

1406:    Level: intermediate

1408:    Notes: A should be obtained by discretizing the vector valued Poisson problem with linear finite elements

1410: .seealso:
1411: @*/
1412: PetscErrorCode PCHYPRESetAlphaPoissonMatrix(PC pc, Mat A)
1413: {

1420:   PetscTryMethod(pc,"PCHYPRESetPoissonMatrix_C",(PC,Mat,PetscBool),(pc,A,PETSC_TRUE));
1421:   return(0);
1422: }

1424: /*@
1425:  PCHYPRESetBetaPoissonMatrix - Set Poisson matrix

1427:    Collective on PC

1429:    Input Parameters:
1430: +  pc - the preconditioning context
1431: -  A - the matrix

1433:    Level: intermediate

1435:    Notes: A should be obtained by discretizing the Poisson problem with linear finite elements.
1436:           Following HYPRE convention, the scalar Poisson solver of AMS can be turned off by passing NULL.

1438: .seealso:
1439: @*/
1440: PetscErrorCode PCHYPRESetBetaPoissonMatrix(PC pc, Mat A)
1441: {

1446:   if (A) {
1449:   }
1450:   PetscTryMethod(pc,"PCHYPRESetPoissonMatrix_C",(PC,Mat,PetscBool),(pc,A,PETSC_FALSE));
1451:   return(0);
1452: }

1454: static PetscErrorCode PCHYPRESetEdgeConstantVectors_HYPRE(PC pc,Vec ozz, Vec zoz, Vec zzo)
1455: {
1456:   PC_HYPRE           *jac = (PC_HYPRE*)pc->data;
1457:   PetscErrorCode     ierr;

1460:   /* throw away any vector if already set */
1461:   if (jac->constants[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[0]));
1462:   if (jac->constants[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[1]));
1463:   if (jac->constants[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[2]));
1464:   jac->constants[0] = NULL;
1465:   jac->constants[1] = NULL;
1466:   jac->constants[2] = NULL;
1467:   VecHYPRE_IJVectorCreate(ozz,&jac->constants[0]);
1468:   VecHYPRE_IJVectorCopy(ozz,jac->constants[0]);
1469:   VecHYPRE_IJVectorCreate(zoz,&jac->constants[1]);
1470:   VecHYPRE_IJVectorCopy(zoz,jac->constants[1]);
1471:   jac->dim = 2;
1472:   if (zzo) {
1473:     VecHYPRE_IJVectorCreate(zzo,&jac->constants[2]);
1474:     VecHYPRE_IJVectorCopy(zzo,jac->constants[2]);
1475:     jac->dim++;
1476:   }
1477:   return(0);
1478: }

1480: /*@
1481:  PCHYPRESetEdgeConstantVectors - Set the representation of the constant vector fields in edge element basis

1483:    Collective on PC

1485:    Input Parameters:
1486: +  pc - the preconditioning context
1487: -  ozz - vector representing (1,0,0) (or (1,0) in 2D)
1488: -  zoz - vector representing (0,1,0) (or (0,1) in 2D)
1489: -  zzo - vector representing (0,0,1) (use NULL in 2D)

1491:    Level: intermediate

1493:    Notes:

1495: .seealso:
1496: @*/
1497: PetscErrorCode PCHYPRESetEdgeConstantVectors(PC pc, Vec ozz, Vec zoz, Vec zzo)
1498: {

1509:   PetscTryMethod(pc,"PCHYPRESetEdgeConstantVectors_C",(PC,Vec,Vec,Vec),(pc,ozz,zoz,zzo));
1510:   return(0);
1511: }

1513: static PetscErrorCode PCSetCoordinates_HYPRE(PC pc, PetscInt dim, PetscInt nloc, PetscReal *coords)
1514: {
1515:   PC_HYPRE        *jac = (PC_HYPRE*)pc->data;
1516:   Vec             tv;
1517:   PetscInt        i;
1518:   PetscErrorCode  ierr;

1521:   /* throw away any coordinate vector if already set */
1522:   if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[0]));
1523:   if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[1]));
1524:   if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[2]));
1525:   jac->dim = dim;

1527:   /* compute IJ vector for coordinates */
1528:   VecCreate(PetscObjectComm((PetscObject)pc),&tv);
1529:   VecSetType(tv,VECSTANDARD);
1530:   VecSetSizes(tv,nloc,PETSC_DECIDE);
1531:   for (i=0;i<dim;i++) {
1532:     PetscScalar *array;
1533:     PetscInt    j;

1535:     VecHYPRE_IJVectorCreate(tv,&jac->coords[i]);
1536:     VecGetArray(tv,&array);
1537:     for (j=0;j<nloc;j++) {
1538:       array[j] = coords[j*dim+i];
1539:     }
1540:     PetscStackCallStandard(HYPRE_IJVectorSetValues,(jac->coords[i],nloc,NULL,array));
1541:     PetscStackCallStandard(HYPRE_IJVectorAssemble,(jac->coords[i]));
1542:     VecRestoreArray(tv,&array);
1543:   }
1544:   VecDestroy(&tv);
1545:   return(0);
1546: }

1548: /* ---------------------------------------------------------------------------------*/

1550: static PetscErrorCode  PCHYPREGetType_HYPRE(PC pc,const char *name[])
1551: {
1552:   PC_HYPRE *jac = (PC_HYPRE*)pc->data;

1555:   *name = jac->hypre_type;
1556:   return(0);
1557: }

1559: static PetscErrorCode  PCHYPRESetType_HYPRE(PC pc,const char name[])
1560: {
1561:   PC_HYPRE       *jac = (PC_HYPRE*)pc->data;
1563:   PetscBool      flag;

1566:   if (jac->hypre_type) {
1567:     PetscStrcmp(jac->hypre_type,name,&flag);
1568:     if (!flag) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ORDER,"Cannot reset the HYPRE preconditioner type once it has been set");
1569:     return(0);
1570:   } else {
1571:     PetscStrallocpy(name, &jac->hypre_type);
1572:   }

1574:   jac->maxiter         = PETSC_DEFAULT;
1575:   jac->tol             = PETSC_DEFAULT;
1576:   jac->printstatistics = PetscLogPrintInfo;

1578:   PetscStrcmp("pilut",jac->hypre_type,&flag);
1579:   if (flag) {
1580:     PetscStackCallStandard(HYPRE_ParCSRPilutCreate,(jac->comm_hypre,&jac->hsolver));
1581:     pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Pilut;
1582:     pc->ops->view           = PCView_HYPRE_Pilut;
1583:     jac->destroy            = HYPRE_ParCSRPilutDestroy;
1584:     jac->setup              = HYPRE_ParCSRPilutSetup;
1585:     jac->solve              = HYPRE_ParCSRPilutSolve;
1586:     jac->factorrowsize      = PETSC_DEFAULT;
1587:     return(0);
1588:   }
1589:   PetscStrcmp("parasails",jac->hypre_type,&flag);
1590:   if (flag) {
1591:     PetscStackCallStandard(HYPRE_ParaSailsCreate,(jac->comm_hypre,&jac->hsolver));
1592:     pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ParaSails;
1593:     pc->ops->view           = PCView_HYPRE_ParaSails;
1594:     jac->destroy            = HYPRE_ParaSailsDestroy;
1595:     jac->setup              = HYPRE_ParaSailsSetup;
1596:     jac->solve              = HYPRE_ParaSailsSolve;
1597:     /* initialize */
1598:     jac->nlevels    = 1;
1599:     jac->threshhold = .1;
1600:     jac->filter     = .1;
1601:     jac->loadbal    = 0;
1602:     if (PetscLogPrintInfo) jac->logging = (int) PETSC_TRUE;
1603:     else jac->logging = (int) PETSC_FALSE;

1605:     jac->ruse = (int) PETSC_FALSE;
1606:     jac->symt = 0;
1607:     PetscStackCallStandard(HYPRE_ParaSailsSetParams,(jac->hsolver,jac->threshhold,jac->nlevels));
1608:     PetscStackCallStandard(HYPRE_ParaSailsSetFilter,(jac->hsolver,jac->filter));
1609:     PetscStackCallStandard(HYPRE_ParaSailsSetLoadbal,(jac->hsolver,jac->loadbal));
1610:     PetscStackCallStandard(HYPRE_ParaSailsSetLogging,(jac->hsolver,jac->logging));
1611:     PetscStackCallStandard(HYPRE_ParaSailsSetReuse,(jac->hsolver,jac->ruse));
1612:     PetscStackCallStandard(HYPRE_ParaSailsSetSym,(jac->hsolver,jac->symt));
1613:     return(0);
1614:   }
1615:   PetscStrcmp("boomeramg",jac->hypre_type,&flag);
1616:   if (flag) {
1617:     HYPRE_BoomerAMGCreate(&jac->hsolver);
1618:     pc->ops->setfromoptions  = PCSetFromOptions_HYPRE_BoomerAMG;
1619:     pc->ops->view            = PCView_HYPRE_BoomerAMG;
1620:     pc->ops->applytranspose  = PCApplyTranspose_HYPRE_BoomerAMG;
1621:     pc->ops->applyrichardson = PCApplyRichardson_HYPRE_BoomerAMG;
1622:     jac->destroy             = HYPRE_BoomerAMGDestroy;
1623:     jac->setup               = HYPRE_BoomerAMGSetup;
1624:     jac->solve               = HYPRE_BoomerAMGSolve;
1625:     jac->applyrichardson     = PETSC_FALSE;
1626:     /* these defaults match the hypre defaults */
1627:     jac->cycletype        = 1;
1628:     jac->maxlevels        = 25;
1629:     jac->maxiter          = 1;
1630:     jac->tol              = 0.0; /* tolerance of zero indicates use as preconditioner (suppresses convergence errors) */
1631:     jac->truncfactor      = 0.0;
1632:     jac->strongthreshold  = .25;
1633:     jac->maxrowsum        = .9;
1634:     jac->coarsentype      = 6;
1635:     jac->measuretype      = 0;
1636:     jac->gridsweeps[0]    = jac->gridsweeps[1] = jac->gridsweeps[2] = 1;
1637:     jac->smoothtype       = -1; /* Not set by default */
1638:     jac->smoothnumlevels  = 25;
1639:     jac->eu_level         = 0;
1640:     jac->eu_droptolerance = 0;
1641:     jac->eu_bj            = 0;
1642:     jac->relaxtype[0]     = jac->relaxtype[1] = 6; /* Defaults to SYMMETRIC since in PETSc we are using a a PC - most likely with CG */
1643:     jac->relaxtype[2]     = 9; /*G.E. */
1644:     jac->relaxweight      = 1.0;
1645:     jac->outerrelaxweight = 1.0;
1646:     jac->relaxorder       = 1;
1647:     jac->interptype       = 0;
1648:     jac->agg_nl           = 0;
1649:     jac->pmax             = 0;
1650:     jac->truncfactor      = 0.0;
1651:     jac->agg_num_paths    = 1;

1653:     jac->nodal_coarsen      = 0;
1654:     jac->nodal_relax        = PETSC_FALSE;
1655:     jac->nodal_relax_levels = 1;
1656:     PetscStackCallStandard(HYPRE_BoomerAMGSetCycleType,(jac->hsolver,jac->cycletype));
1657:     PetscStackCallStandard(HYPRE_BoomerAMGSetMaxLevels,(jac->hsolver,jac->maxlevels));
1658:     PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,jac->maxiter));
1659:     PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,jac->tol));
1660:     PetscStackCallStandard(HYPRE_BoomerAMGSetTruncFactor,(jac->hsolver,jac->truncfactor));
1661:     PetscStackCallStandard(HYPRE_BoomerAMGSetStrongThreshold,(jac->hsolver,jac->strongthreshold));
1662:     PetscStackCallStandard(HYPRE_BoomerAMGSetMaxRowSum,(jac->hsolver,jac->maxrowsum));
1663:     PetscStackCallStandard(HYPRE_BoomerAMGSetCoarsenType,(jac->hsolver,jac->coarsentype));
1664:     PetscStackCallStandard(HYPRE_BoomerAMGSetMeasureType,(jac->hsolver,jac->measuretype));
1665:     PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxOrder,(jac->hsolver, jac->relaxorder));
1666:     PetscStackCallStandard(HYPRE_BoomerAMGSetInterpType,(jac->hsolver,jac->interptype));
1667:     PetscStackCallStandard(HYPRE_BoomerAMGSetAggNumLevels,(jac->hsolver,jac->agg_nl));
1668:     PetscStackCallStandard(HYPRE_BoomerAMGSetPMaxElmts,(jac->hsolver,jac->pmax));
1669:     PetscStackCallStandard(HYPRE_BoomerAMGSetNumPaths,(jac->hsolver,jac->agg_num_paths));
1670:     PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxType,(jac->hsolver, jac->relaxtype[0]));  /*defaults coarse to 9*/
1671:     PetscStackCallStandard(HYPRE_BoomerAMGSetNumSweeps,(jac->hsolver, jac->gridsweeps[0])); /*defaults coarse to 1 */
1672:     return(0);
1673:   }
1674:   PetscStrcmp("ams",jac->hypre_type,&flag);
1675:   if (flag) {
1676:     HYPRE_AMSCreate(&jac->hsolver);
1677:     pc->ops->setfromoptions  = PCSetFromOptions_HYPRE_AMS;
1678:     pc->ops->view            = PCView_HYPRE_AMS;
1679:     jac->destroy             = HYPRE_AMSDestroy;
1680:     jac->setup               = HYPRE_AMSSetup;
1681:     jac->solve               = HYPRE_AMSSolve;
1682:     jac->coords[0]           = NULL;
1683:     jac->coords[1]           = NULL;
1684:     jac->coords[2]           = NULL;
1685:     /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE AMS */
1686:     jac->as_print           = 0;
1687:     jac->as_max_iter        = 1; /* used as a preconditioner */
1688:     jac->as_tol             = 0.; /* used as a preconditioner */
1689:     jac->ams_cycle_type     = 13;
1690:     /* Smoothing options */
1691:     jac->as_relax_type      = 2;
1692:     jac->as_relax_times     = 1;
1693:     jac->as_relax_weight    = 1.0;
1694:     jac->as_omega           = 1.0;
1695:     /* Vector valued Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
1696:     jac->as_amg_alpha_opts[0] = 10;
1697:     jac->as_amg_alpha_opts[1] = 1;
1698:     jac->as_amg_alpha_opts[2] = 6;
1699:     jac->as_amg_alpha_opts[3] = 6;
1700:     jac->as_amg_alpha_opts[4] = 4;
1701:     jac->as_amg_alpha_theta   = 0.25;
1702:     /* Scalar Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
1703:     jac->as_amg_beta_opts[0] = 10;
1704:     jac->as_amg_beta_opts[1] = 1;
1705:     jac->as_amg_beta_opts[2] = 6;
1706:     jac->as_amg_beta_opts[3] = 6;
1707:     jac->as_amg_beta_opts[4] = 4;
1708:     jac->as_amg_beta_theta   = 0.25;
1709:     PetscStackCallStandard(HYPRE_AMSSetPrintLevel,(jac->hsolver,jac->as_print));
1710:     PetscStackCallStandard(HYPRE_AMSSetMaxIter,(jac->hsolver,jac->as_max_iter));
1711:     PetscStackCallStandard(HYPRE_AMSSetCycleType,(jac->hsolver,jac->ams_cycle_type));
1712:     PetscStackCallStandard(HYPRE_AMSSetTol,(jac->hsolver,jac->as_tol));
1713:     PetscStackCallStandard(HYPRE_AMSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type,
1714:                                                                       jac->as_relax_times,
1715:                                                                       jac->as_relax_weight,
1716:                                                                       jac->as_omega));
1717:     PetscStackCallStandard(HYPRE_AMSSetAlphaAMGOptions,(jac->hsolver,jac->as_amg_alpha_opts[0],       /* AMG coarsen type */
1718:                                                                      jac->as_amg_alpha_opts[1],       /* AMG agg_levels */
1719:                                                                      jac->as_amg_alpha_opts[2],       /* AMG relax_type */
1720:                                                                      jac->as_amg_alpha_theta,
1721:                                                                      jac->as_amg_alpha_opts[3],       /* AMG interp_type */
1722:                                                                      jac->as_amg_alpha_opts[4]));     /* AMG Pmax */
1723:     PetscStackCallStandard(HYPRE_AMSSetBetaAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0],       /* AMG coarsen type */
1724:                                                                     jac->as_amg_beta_opts[1],       /* AMG agg_levels */
1725:                                                                     jac->as_amg_beta_opts[2],       /* AMG relax_type */
1726:                                                                     jac->as_amg_beta_theta,
1727:                                                                     jac->as_amg_beta_opts[3],       /* AMG interp_type */
1728:                                                                     jac->as_amg_beta_opts[4]));     /* AMG Pmax */
1729:     /* Zero conductivity */
1730:     jac->ams_beta_is_zero      = PETSC_FALSE;
1731:     jac->ams_beta_is_zero_part = PETSC_FALSE;
1732:     return(0);
1733:   }
1734:   PetscStrcmp("ads",jac->hypre_type,&flag);
1735:   if (flag) {
1736:     HYPRE_ADSCreate(&jac->hsolver);
1737:     pc->ops->setfromoptions  = PCSetFromOptions_HYPRE_ADS;
1738:     pc->ops->view            = PCView_HYPRE_ADS;
1739:     jac->destroy             = HYPRE_ADSDestroy;
1740:     jac->setup               = HYPRE_ADSSetup;
1741:     jac->solve               = HYPRE_ADSSolve;
1742:     jac->coords[0]           = NULL;
1743:     jac->coords[1]           = NULL;
1744:     jac->coords[2]           = NULL;
1745:     /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE ADS */
1746:     jac->as_print           = 0;
1747:     jac->as_max_iter        = 1; /* used as a preconditioner */
1748:     jac->as_tol             = 0.; /* used as a preconditioner */
1749:     jac->ads_cycle_type     = 13;
1750:     /* Smoothing options */
1751:     jac->as_relax_type      = 2;
1752:     jac->as_relax_times     = 1;
1753:     jac->as_relax_weight    = 1.0;
1754:     jac->as_omega           = 1.0;
1755:     /* AMS solver parameters: cycle_type, coarsen type, agg_levels, relax_type, interp_type, Pmax */
1756:     jac->ams_cycle_type       = 14;
1757:     jac->as_amg_alpha_opts[0] = 10;
1758:     jac->as_amg_alpha_opts[1] = 1;
1759:     jac->as_amg_alpha_opts[2] = 6;
1760:     jac->as_amg_alpha_opts[3] = 6;
1761:     jac->as_amg_alpha_opts[4] = 4;
1762:     jac->as_amg_alpha_theta   = 0.25;
1763:     /* Vector Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
1764:     jac->as_amg_beta_opts[0] = 10;
1765:     jac->as_amg_beta_opts[1] = 1;
1766:     jac->as_amg_beta_opts[2] = 6;
1767:     jac->as_amg_beta_opts[3] = 6;
1768:     jac->as_amg_beta_opts[4] = 4;
1769:     jac->as_amg_beta_theta   = 0.25;
1770:     PetscStackCallStandard(HYPRE_ADSSetPrintLevel,(jac->hsolver,jac->as_print));
1771:     PetscStackCallStandard(HYPRE_ADSSetMaxIter,(jac->hsolver,jac->as_max_iter));
1772:     PetscStackCallStandard(HYPRE_ADSSetCycleType,(jac->hsolver,jac->ams_cycle_type));
1773:     PetscStackCallStandard(HYPRE_ADSSetTol,(jac->hsolver,jac->as_tol));
1774:     PetscStackCallStandard(HYPRE_ADSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type,
1775:                                                                       jac->as_relax_times,
1776:                                                                       jac->as_relax_weight,
1777:                                                                       jac->as_omega));
1778:     PetscStackCallStandard(HYPRE_ADSSetAMSOptions,(jac->hsolver,jac->ams_cycle_type,             /* AMG coarsen type */
1779:                                                                 jac->as_amg_alpha_opts[0],       /* AMG coarsen type */
1780:                                                                 jac->as_amg_alpha_opts[1],       /* AMG agg_levels */
1781:                                                                 jac->as_amg_alpha_opts[2],       /* AMG relax_type */
1782:                                                                 jac->as_amg_alpha_theta,
1783:                                                                 jac->as_amg_alpha_opts[3],       /* AMG interp_type */
1784:                                                                 jac->as_amg_alpha_opts[4]));     /* AMG Pmax */
1785:     PetscStackCallStandard(HYPRE_ADSSetAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0],       /* AMG coarsen type */
1786:                                                                 jac->as_amg_beta_opts[1],       /* AMG agg_levels */
1787:                                                                 jac->as_amg_beta_opts[2],       /* AMG relax_type */
1788:                                                                 jac->as_amg_beta_theta,
1789:                                                                 jac->as_amg_beta_opts[3],       /* AMG interp_type */
1790:                                                                 jac->as_amg_beta_opts[4]));     /* AMG Pmax */
1791:     return(0);
1792:   }
1793:   PetscFree(jac->hypre_type);

1795:   jac->hypre_type = NULL;
1796:   SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown HYPRE preconditioner %s; Choices are pilut, parasails, boomeramg, ams",name);
1797:   return(0);
1798: }

1800: /*
1801:     It only gets here if the HYPRE type has not been set before the call to
1802:    ...SetFromOptions() which actually is most of the time
1803: */
1804: static PetscErrorCode PCSetFromOptions_HYPRE(PetscOptionItems *PetscOptionsObject,PC pc)
1805: {
1807:   PetscInt       indx;
1808:   const char     *type[] = {"pilut","parasails","boomeramg","ams","ads"};
1809:   PetscBool      flg;

1812:   PetscOptionsHead(PetscOptionsObject,"HYPRE preconditioner options");
1813:   PetscOptionsEList("-pc_hypre_type","HYPRE preconditioner type","PCHYPRESetType",type,4,"boomeramg",&indx,&flg);
1814:   if (flg) {
1815:     PCHYPRESetType_HYPRE(pc,type[indx]);
1816:   } else {
1817:     PCHYPRESetType_HYPRE(pc,"boomeramg");
1818:   }
1819:   if (pc->ops->setfromoptions) {
1820:     pc->ops->setfromoptions(PetscOptionsObject,pc);
1821:   }
1822:   PetscOptionsTail();
1823:   return(0);
1824: }

1826: /*@C
1827:      PCHYPRESetType - Sets which hypre preconditioner you wish to use

1829:    Input Parameters:
1830: +     pc - the preconditioner context
1831: -     name - either  pilut, parasails, boomeramg, ams, ads

1833:    Options Database Keys:
1834:    -pc_hypre_type - One of pilut, parasails, boomeramg, ams, ads

1836:    Level: intermediate

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

1841: @*/
1842: PetscErrorCode  PCHYPRESetType(PC pc,const char name[])
1843: {

1849:   PetscTryMethod(pc,"PCHYPRESetType_C",(PC,const char[]),(pc,name));
1850:   return(0);
1851: }

1853: /*@C
1854:      PCHYPREGetType - Gets which hypre preconditioner you are using

1856:    Input Parameter:
1857: .     pc - the preconditioner context

1859:    Output Parameter:
1860: .     name - either  pilut, parasails, boomeramg, ams, ads

1862:    Level: intermediate

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

1867: @*/
1868: PetscErrorCode  PCHYPREGetType(PC pc,const char *name[])
1869: {

1875:   PetscTryMethod(pc,"PCHYPREGetType_C",(PC,const char*[]),(pc,name));
1876:   return(0);
1877: }

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

1882:    Options Database Keys:
1883: +   -pc_hypre_type - One of pilut, parasails, boomeramg, ams, ads
1884: -   Too many others to list, run with -pc_type hypre -pc_hypre_type XXX -help to see options for the XXX
1885:           preconditioner

1887:    Level: intermediate

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

1893:           The options -pc_hypre_boomeramg_max_iter and -pc_hypre_boomeramg_tol refer to the number of iterations
1894:           (V-cycles) and tolerance that boomeramg does EACH time it is called. So for example, if
1895:           -pc_hypre_boomeramg_max_iter is set to 2 then 2-V-cycles are being used to define the preconditioner
1896:           (-pc_hypre_boomeramg_tol should be set to 0.0 - the default - to strictly use a fixed number of
1897:           iterations per hypre call). -ksp_max_it and -ksp_rtol STILL determine the total number of iterations
1898:           and tolerance for the Krylov solver. For example, if -pc_hypre_boomeramg_max_iter is 2 and -ksp_max_it is 10
1899:           then AT MOST twenty V-cycles of boomeramg will be called.

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

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

1911:           MatSetNearNullSpace() - if you provide a near null space to your matrix it is ignored by hypre UNLESS you also use
1912:           the two options:
1913: +   -pc_hypre_boomeramg_nodal_coarsen <n> - where n is from 1 to 6 (see HYPRE_BOOMERAMGSetNodal())
1914: -   -pc_hypre_boomeramg_vec_interp_variant <v> where v is from 1 to 3 (see HYPRE_BoomerAMGSetInterpVecVariant())

1916:           Depending on the linear system you may see the same or different convergence depending on the values you use.

1918:           See PCPFMG for access to the hypre Struct PFMG solver

1920: .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PC,
1921:            PCHYPRESetType(), PCPFMG

1923: M*/

1925: PETSC_EXTERN PetscErrorCode PCCreate_HYPRE(PC pc)
1926: {
1927:   PC_HYPRE       *jac;

1931:   PetscNewLog(pc,&jac);

1933:   pc->data                = jac;
1934:   pc->ops->reset          = PCReset_HYPRE;
1935:   pc->ops->destroy        = PCDestroy_HYPRE;
1936:   pc->ops->setfromoptions = PCSetFromOptions_HYPRE;
1937:   pc->ops->setup          = PCSetUp_HYPRE;
1938:   pc->ops->apply          = PCApply_HYPRE;
1939:   jac->comm_hypre         = MPI_COMM_NULL;
1940:   /* duplicate communicator for hypre */
1941:   MPI_Comm_dup(PetscObjectComm((PetscObject)pc),&(jac->comm_hypre));
1942:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetType_C",PCHYPRESetType_HYPRE);
1943:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPREGetType_C",PCHYPREGetType_HYPRE);
1944:   PetscObjectComposeFunction((PetscObject)pc,"PCSetCoordinates_C",PCSetCoordinates_HYPRE);
1945:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteGradient_C",PCHYPRESetDiscreteGradient_HYPRE);
1946:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteCurl_C",PCHYPRESetDiscreteCurl_HYPRE);
1947:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetInterpolations_C",PCHYPRESetInterpolations_HYPRE);
1948:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetEdgeConstantVectors_C",PCHYPRESetEdgeConstantVectors_HYPRE);
1949:   PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetPoissonMatrix_C",PCHYPRESetPoissonMatrix_HYPRE);
1950:   return(0);
1951: }

1953: /* ---------------------------------------------------------------------------------------------------------------------------------*/

1955: typedef struct {
1956:   MPI_Comm           hcomm;        /* does not share comm with HYPRE_StructMatrix because need to create solver before getting matrix */
1957:   HYPRE_StructSolver hsolver;

1959:   /* keep copy of PFMG options used so may view them */
1960:   PetscInt its;
1961:   double   tol;
1962:   PetscInt relax_type;
1963:   PetscInt rap_type;
1964:   PetscInt num_pre_relax,num_post_relax;
1965:   PetscInt max_levels;
1966: } PC_PFMG;

1968: PetscErrorCode PCDestroy_PFMG(PC pc)
1969: {
1971:   PC_PFMG        *ex = (PC_PFMG*) pc->data;

1974:   if (ex->hsolver) PetscStackCallStandard(HYPRE_StructPFMGDestroy,(ex->hsolver));
1975:   MPI_Comm_free(&ex->hcomm);
1976:   PetscFree(pc->data);
1977:   return(0);
1978: }

1980: static const char *PFMGRelaxType[] = {"Jacobi","Weighted-Jacobi","symmetric-Red/Black-Gauss-Seidel","Red/Black-Gauss-Seidel"};
1981: static const char *PFMGRAPType[] = {"Galerkin","non-Galerkin"};

1983: PetscErrorCode PCView_PFMG(PC pc,PetscViewer viewer)
1984: {
1986:   PetscBool      iascii;
1987:   PC_PFMG        *ex = (PC_PFMG*) pc->data;

1990:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
1991:   if (iascii) {
1992:     PetscViewerASCIIPrintf(viewer,"  HYPRE PFMG preconditioning\n");
1993:     PetscViewerASCIIPrintf(viewer,"    max iterations %d\n",ex->its);
1994:     PetscViewerASCIIPrintf(viewer,"    tolerance %g\n",ex->tol);
1995:     PetscViewerASCIIPrintf(viewer,"    relax type %s\n",PFMGRelaxType[ex->relax_type]);
1996:     PetscViewerASCIIPrintf(viewer,"    RAP type %s\n",PFMGRAPType[ex->rap_type]);
1997:     PetscViewerASCIIPrintf(viewer,"    number pre-relax %d post-relax %d\n",ex->num_pre_relax,ex->num_post_relax);
1998:     PetscViewerASCIIPrintf(viewer,"    max levels %d\n",ex->max_levels);
1999:   }
2000:   return(0);
2001: }

2003: PetscErrorCode PCSetFromOptions_PFMG(PetscOptionItems *PetscOptionsObject,PC pc)
2004: {
2006:   PC_PFMG        *ex = (PC_PFMG*) pc->data;
2007:   PetscBool      flg = PETSC_FALSE;

2010:   PetscOptionsHead(PetscOptionsObject,"PFMG options");
2011:   PetscOptionsBool("-pc_pfmg_print_statistics","Print statistics","HYPRE_StructPFMGSetPrintLevel",flg,&flg,NULL);
2012:   if (flg) {
2013:     int level=3;
2014:     PetscStackCallStandard(HYPRE_StructPFMGSetPrintLevel,(ex->hsolver,level));
2015:   }
2016:   PetscOptionsInt("-pc_pfmg_its","Number of iterations of PFMG to use as preconditioner","HYPRE_StructPFMGSetMaxIter",ex->its,&ex->its,NULL);
2017:   PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,(ex->hsolver,ex->its));
2018:   PetscOptionsInt("-pc_pfmg_num_pre_relax","Number of smoothing steps before coarse grid","HYPRE_StructPFMGSetNumPreRelax",ex->num_pre_relax,&ex->num_pre_relax,NULL);
2019:   PetscStackCallStandard(HYPRE_StructPFMGSetNumPreRelax,(ex->hsolver,ex->num_pre_relax));
2020:   PetscOptionsInt("-pc_pfmg_num_post_relax","Number of smoothing steps after coarse grid","HYPRE_StructPFMGSetNumPostRelax",ex->num_post_relax,&ex->num_post_relax,NULL);
2021:   PetscStackCallStandard(HYPRE_StructPFMGSetNumPostRelax,(ex->hsolver,ex->num_post_relax));

2023:   PetscOptionsInt("-pc_pfmg_max_levels","Max Levels for MG hierarchy","HYPRE_StructPFMGSetMaxLevels",ex->max_levels,&ex->max_levels,NULL);
2024:   PetscStackCallStandard(HYPRE_StructPFMGSetMaxLevels,(ex->hsolver,ex->max_levels));

2026:   PetscOptionsReal("-pc_pfmg_tol","Tolerance of PFMG","HYPRE_StructPFMGSetTol",ex->tol,&ex->tol,NULL);
2027:   PetscStackCallStandard(HYPRE_StructPFMGSetTol,(ex->hsolver,ex->tol));
2028:   PetscOptionsEList("-pc_pfmg_relax_type","Relax type for the up and down cycles","HYPRE_StructPFMGSetRelaxType",PFMGRelaxType,ALEN(PFMGRelaxType),PFMGRelaxType[ex->relax_type],&ex->relax_type,NULL);
2029:   PetscStackCallStandard(HYPRE_StructPFMGSetRelaxType,(ex->hsolver, ex->relax_type));
2030:   PetscOptionsEList("-pc_pfmg_rap_type","RAP type","HYPRE_StructPFMGSetRAPType",PFMGRAPType,ALEN(PFMGRAPType),PFMGRAPType[ex->rap_type],&ex->rap_type,NULL);
2031:   PetscStackCallStandard(HYPRE_StructPFMGSetRAPType,(ex->hsolver, ex->rap_type));
2032:   PetscOptionsTail();
2033:   return(0);
2034: }

2036: PetscErrorCode PCApply_PFMG(PC pc,Vec x,Vec y)
2037: {
2038:   PetscErrorCode    ierr;
2039:   PC_PFMG           *ex = (PC_PFMG*) pc->data;
2040:   PetscScalar       *yy;
2041:   const PetscScalar *xx;
2042:   PetscInt          ilower[3],iupper[3];
2043:   Mat_HYPREStruct   *mx = (Mat_HYPREStruct*)(pc->pmat->data);

2046:   PetscCitationsRegister(hypreCitation,&cite);
2047:   DMDAGetCorners(mx->da,&ilower[0],&ilower[1],&ilower[2],&iupper[0],&iupper[1],&iupper[2]);
2048:   iupper[0] += ilower[0] - 1;
2049:   iupper[1] += ilower[1] - 1;
2050:   iupper[2] += ilower[2] - 1;

2052:   /* copy x values over to hypre */
2053:   PetscStackCallStandard(HYPRE_StructVectorSetConstantValues,(mx->hb,0.0));
2054:   VecGetArrayRead(x,&xx);
2055:   PetscStackCallStandard(HYPRE_StructVectorSetBoxValues,(mx->hb,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,(PetscScalar*)xx));
2056:   VecRestoreArrayRead(x,&xx);
2057:   PetscStackCallStandard(HYPRE_StructVectorAssemble,(mx->hb));
2058:   PetscStackCallStandard(HYPRE_StructPFMGSolve,(ex->hsolver,mx->hmat,mx->hb,mx->hx));

2060:   /* copy solution values back to PETSc */
2061:   VecGetArray(y,&yy);
2062:   PetscStackCallStandard(HYPRE_StructVectorGetBoxValues,(mx->hx,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,yy));
2063:   VecRestoreArray(y,&yy);
2064:   return(0);
2065: }

2067: static PetscErrorCode PCApplyRichardson_PFMG(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its,PetscBool guesszero,PetscInt *outits,PCRichardsonConvergedReason *reason)
2068: {
2069:   PC_PFMG        *jac = (PC_PFMG*)pc->data;
2071:   PetscInt       oits;

2074:   PetscCitationsRegister(hypreCitation,&cite);
2075:   PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,(jac->hsolver,its*jac->its));
2076:   PetscStackCallStandard(HYPRE_StructPFMGSetTol,(jac->hsolver,rtol));

2078:   PCApply_PFMG(pc,b,y);
2079:   PetscStackCallStandard(HYPRE_StructPFMGGetNumIterations,(jac->hsolver,(HYPRE_Int *)&oits));
2080:   *outits = oits;
2081:   if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
2082:   else             *reason = PCRICHARDSON_CONVERGED_RTOL;
2083:   PetscStackCallStandard(HYPRE_StructPFMGSetTol,(jac->hsolver,jac->tol));
2084:   PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,(jac->hsolver,jac->its));
2085:   return(0);
2086: }


2089: PetscErrorCode PCSetUp_PFMG(PC pc)
2090: {
2091:   PetscErrorCode  ierr;
2092:   PC_PFMG         *ex = (PC_PFMG*) pc->data;
2093:   Mat_HYPREStruct *mx = (Mat_HYPREStruct*)(pc->pmat->data);
2094:   PetscBool       flg;

2097:   PetscObjectTypeCompare((PetscObject)pc->pmat,MATHYPRESTRUCT,&flg);
2098:   if (!flg) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Must use MATHYPRESTRUCT with this preconditioner");

2100:   /* create the hypre solver object and set its information */
2101:   if (ex->hsolver) PetscStackCallStandard(HYPRE_StructPFMGDestroy,(ex->hsolver));
2102:   PetscStackCallStandard(HYPRE_StructPFMGCreate,(ex->hcomm,&ex->hsolver));
2103:   PetscStackCallStandard(HYPRE_StructPFMGSetup,(ex->hsolver,mx->hmat,mx->hb,mx->hx));
2104:   PetscStackCallStandard(HYPRE_StructPFMGSetZeroGuess,(ex->hsolver));
2105:   return(0);
2106: }

2108: /*MC
2109:      PCPFMG - the hypre PFMG multigrid solver

2111:    Level: advanced

2113:    Options Database:
2114: + -pc_pfmg_its <its> number of iterations of PFMG to use as preconditioner
2115: . -pc_pfmg_num_pre_relax <steps> number of smoothing steps before coarse grid
2116: . -pc_pfmg_num_post_relax <steps> number of smoothing steps after coarse grid
2117: . -pc_pfmg_tol <tol> tolerance of PFMG
2118: . -pc_pfmg_relax_type -relaxation type for the up and down cycles, one of Jacobi,Weighted-Jacobi,symmetric-Red/Black-Gauss-Seidel,Red/Black-Gauss-Seidel
2119: - -pc_pfmg_rap_type - type of coarse matrix generation, one of Galerkin,non-Galerkin

2121:    Notes:  This is for CELL-centered descretizations

2123:            This must be used with the MATHYPRESTRUCT matrix type.
2124:            This is less general than in hypre, it supports only one block per process defined by a PETSc DMDA.

2126: .seealso:  PCMG, MATHYPRESTRUCT
2127: M*/

2129: PETSC_EXTERN PetscErrorCode PCCreate_PFMG(PC pc)
2130: {
2132:   PC_PFMG        *ex;

2135:   PetscNew(&ex); \
2136:   pc->data = ex;

2138:   ex->its            = 1;
2139:   ex->tol            = 1.e-8;
2140:   ex->relax_type     = 1;
2141:   ex->rap_type       = 0;
2142:   ex->num_pre_relax  = 1;
2143:   ex->num_post_relax = 1;
2144:   ex->max_levels     = 0;

2146:   pc->ops->setfromoptions  = PCSetFromOptions_PFMG;
2147:   pc->ops->view            = PCView_PFMG;
2148:   pc->ops->destroy         = PCDestroy_PFMG;
2149:   pc->ops->apply           = PCApply_PFMG;
2150:   pc->ops->applyrichardson = PCApplyRichardson_PFMG;
2151:   pc->ops->setup           = PCSetUp_PFMG;

2153:   MPI_Comm_dup(PetscObjectComm((PetscObject)pc),&(ex->hcomm));
2154:   PetscStackCallStandard(HYPRE_StructPFMGCreate,(ex->hcomm,&ex->hsolver));
2155:   return(0);
2156: }

2158: /* ---------------------------------------------------------------------------------------------------------------------------------------------------*/

2160: /* we know we are working with a HYPRE_SStructMatrix */
2161: typedef struct {
2162:   MPI_Comm            hcomm;       /* does not share comm with HYPRE_SStructMatrix because need to create solver before getting matrix */
2163:   HYPRE_SStructSolver ss_solver;

2165:   /* keep copy of SYSPFMG options used so may view them */
2166:   PetscInt its;
2167:   double   tol;
2168:   PetscInt relax_type;
2169:   PetscInt num_pre_relax,num_post_relax;
2170: } PC_SysPFMG;

2172: PetscErrorCode PCDestroy_SysPFMG(PC pc)
2173: {
2175:   PC_SysPFMG     *ex = (PC_SysPFMG*) pc->data;

2178:   if (ex->ss_solver) PetscStackCallStandard(HYPRE_SStructSysPFMGDestroy,(ex->ss_solver));
2179:   MPI_Comm_free(&ex->hcomm);
2180:   PetscFree(pc->data);
2181:   return(0);
2182: }

2184: static const char *SysPFMGRelaxType[] = {"Weighted-Jacobi","Red/Black-Gauss-Seidel"};

2186: PetscErrorCode PCView_SysPFMG(PC pc,PetscViewer viewer)
2187: {
2189:   PetscBool      iascii;
2190:   PC_SysPFMG     *ex = (PC_SysPFMG*) pc->data;

2193:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
2194:   if (iascii) {
2195:     PetscViewerASCIIPrintf(viewer,"  HYPRE SysPFMG preconditioning\n");
2196:     PetscViewerASCIIPrintf(viewer,"  max iterations %d\n",ex->its);
2197:     PetscViewerASCIIPrintf(viewer,"  tolerance %g\n",ex->tol);
2198:     PetscViewerASCIIPrintf(viewer,"  relax type %s\n",PFMGRelaxType[ex->relax_type]);
2199:     PetscViewerASCIIPrintf(viewer,"  number pre-relax %d post-relax %d\n",ex->num_pre_relax,ex->num_post_relax);
2200:   }
2201:   return(0);
2202: }

2204: PetscErrorCode PCSetFromOptions_SysPFMG(PetscOptionItems *PetscOptionsObject,PC pc)
2205: {
2207:   PC_SysPFMG     *ex = (PC_SysPFMG*) pc->data;
2208:   PetscBool      flg = PETSC_FALSE;

2211:   PetscOptionsHead(PetscOptionsObject,"SysPFMG options");
2212:   PetscOptionsBool("-pc_syspfmg_print_statistics","Print statistics","HYPRE_SStructSysPFMGSetPrintLevel",flg,&flg,NULL);
2213:   if (flg) {
2214:     int level=3;
2215:     PetscStackCallStandard(HYPRE_SStructSysPFMGSetPrintLevel,(ex->ss_solver,level));
2216:   }
2217:   PetscOptionsInt("-pc_syspfmg_its","Number of iterations of SysPFMG to use as preconditioner","HYPRE_SStructSysPFMGSetMaxIter",ex->its,&ex->its,NULL);
2218:   PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,(ex->ss_solver,ex->its));
2219:   PetscOptionsInt("-pc_syspfmg_num_pre_relax","Number of smoothing steps before coarse grid","HYPRE_SStructSysPFMGSetNumPreRelax",ex->num_pre_relax,&ex->num_pre_relax,NULL);
2220:   PetscStackCallStandard(HYPRE_SStructSysPFMGSetNumPreRelax,(ex->ss_solver,ex->num_pre_relax));
2221:   PetscOptionsInt("-pc_syspfmg_num_post_relax","Number of smoothing steps after coarse grid","HYPRE_SStructSysPFMGSetNumPostRelax",ex->num_post_relax,&ex->num_post_relax,NULL);
2222:   PetscStackCallStandard(HYPRE_SStructSysPFMGSetNumPostRelax,(ex->ss_solver,ex->num_post_relax));

2224:   PetscOptionsReal("-pc_syspfmg_tol","Tolerance of SysPFMG","HYPRE_SStructSysPFMGSetTol",ex->tol,&ex->tol,NULL);
2225:   PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,(ex->ss_solver,ex->tol));
2226:   PetscOptionsEList("-pc_syspfmg_relax_type","Relax type for the up and down cycles","HYPRE_SStructSysPFMGSetRelaxType",SysPFMGRelaxType,4,SysPFMGRelaxType[ex->relax_type],&ex->relax_type,NULL);
2227:   PetscStackCallStandard(HYPRE_SStructSysPFMGSetRelaxType,(ex->ss_solver, ex->relax_type));
2228:   PetscOptionsTail();
2229:   return(0);
2230: }

2232: PetscErrorCode PCApply_SysPFMG(PC pc,Vec x,Vec y)
2233: {
2234:   PetscErrorCode    ierr;
2235:   PC_SysPFMG        *ex = (PC_SysPFMG*) pc->data;
2236:   PetscScalar       *yy;
2237:   const PetscScalar *xx;
2238:   PetscInt          ilower[3],iupper[3];
2239:   Mat_HYPRESStruct  *mx     = (Mat_HYPRESStruct*)(pc->pmat->data);
2240:   PetscInt          ordering= mx->dofs_order;
2241:   PetscInt          nvars   = mx->nvars;
2242:   PetscInt          part    = 0;
2243:   PetscInt          size;
2244:   PetscInt          i;

2247:   PetscCitationsRegister(hypreCitation,&cite);
2248:   DMDAGetCorners(mx->da,&ilower[0],&ilower[1],&ilower[2],&iupper[0],&iupper[1],&iupper[2]);
2249:   iupper[0] += ilower[0] - 1;
2250:   iupper[1] += ilower[1] - 1;
2251:   iupper[2] += ilower[2] - 1;

2253:   size = 1;
2254:   for (i= 0; i< 3; i++) size *= (iupper[i]-ilower[i]+1);

2256:   /* copy x values over to hypre for variable ordering */
2257:   if (ordering) {
2258:     PetscStackCallStandard(HYPRE_SStructVectorSetConstantValues,(mx->ss_b,0.0));
2259:     VecGetArrayRead(x,&xx);
2260:     for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorSetBoxValues,(mx->ss_b,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,(PetscScalar*)xx+(size*i)));
2261:     VecRestoreArrayRead(x,&xx);
2262:     PetscStackCallStandard(HYPRE_SStructVectorAssemble,(mx->ss_b));
2263:     PetscStackCallStandard(HYPRE_SStructMatrixMatvec,(1.0,mx->ss_mat,mx->ss_b,0.0,mx->ss_x));
2264:     PetscStackCallStandard(HYPRE_SStructSysPFMGSolve,(ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x));

2266:     /* copy solution values back to PETSc */
2267:     VecGetArray(y,&yy);
2268:     for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorGetBoxValues,(mx->ss_x,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,yy+(size*i)));
2269:     VecRestoreArray(y,&yy);
2270:   } else {      /* nodal ordering must be mapped to variable ordering for sys_pfmg */
2271:     PetscScalar *z;
2272:     PetscInt    j, k;

2274:     PetscMalloc1(nvars*size,&z);
2275:     PetscStackCallStandard(HYPRE_SStructVectorSetConstantValues,(mx->ss_b,0.0));
2276:     VecGetArrayRead(x,&xx);

2278:     /* transform nodal to hypre's variable ordering for sys_pfmg */
2279:     for (i= 0; i< size; i++) {
2280:       k= i*nvars;
2281:       for (j= 0; j< nvars; j++) z[j*size+i]= xx[k+j];
2282:     }
2283:     for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorSetBoxValues,(mx->ss_b,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,z+(size*i)));
2284:     VecRestoreArrayRead(x,&xx);
2285:     PetscStackCallStandard(HYPRE_SStructVectorAssemble,(mx->ss_b));
2286:     PetscStackCallStandard(HYPRE_SStructSysPFMGSolve,(ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x));

2288:     /* copy solution values back to PETSc */
2289:     VecGetArray(y,&yy);
2290:     for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorGetBoxValues,(mx->ss_x,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,z+(size*i)));
2291:     /* transform hypre's variable ordering for sys_pfmg to nodal ordering */
2292:     for (i= 0; i< size; i++) {
2293:       k= i*nvars;
2294:       for (j= 0; j< nvars; j++) yy[k+j]= z[j*size+i];
2295:     }
2296:     VecRestoreArray(y,&yy);
2297:     PetscFree(z);
2298:   }
2299:   return(0);
2300: }

2302: static PetscErrorCode PCApplyRichardson_SysPFMG(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its,PetscBool guesszero,PetscInt *outits,PCRichardsonConvergedReason *reason)
2303: {
2304:   PC_SysPFMG     *jac = (PC_SysPFMG*)pc->data;
2306:   PetscInt       oits;

2309:   PetscCitationsRegister(hypreCitation,&cite);
2310:   PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,(jac->ss_solver,its*jac->its));
2311:   PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,(jac->ss_solver,rtol));
2312:   PCApply_SysPFMG(pc,b,y);
2313:   PetscStackCallStandard(HYPRE_SStructSysPFMGGetNumIterations,(jac->ss_solver,(HYPRE_Int *)&oits));
2314:   *outits = oits;
2315:   if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
2316:   else             *reason = PCRICHARDSON_CONVERGED_RTOL;
2317:   PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,(jac->ss_solver,jac->tol));
2318:   PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,(jac->ss_solver,jac->its));
2319:   return(0);
2320: }

2322: PetscErrorCode PCSetUp_SysPFMG(PC pc)
2323: {
2324:   PetscErrorCode   ierr;
2325:   PC_SysPFMG       *ex = (PC_SysPFMG*) pc->data;
2326:   Mat_HYPRESStruct *mx = (Mat_HYPRESStruct*)(pc->pmat->data);
2327:   PetscBool        flg;

2330:   PetscObjectTypeCompare((PetscObject)pc->pmat,MATHYPRESSTRUCT,&flg);
2331:   if (!flg) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Must use MATHYPRESSTRUCT with this preconditioner");

2333:   /* create the hypre sstruct solver object and set its information */
2334:   if (ex->ss_solver) PetscStackCallStandard(HYPRE_SStructSysPFMGDestroy,(ex->ss_solver));
2335:   PetscStackCallStandard(HYPRE_SStructSysPFMGCreate,(ex->hcomm,&ex->ss_solver));
2336:   PetscStackCallStandard(HYPRE_SStructSysPFMGSetZeroGuess,(ex->ss_solver));
2337:   PetscStackCallStandard(HYPRE_SStructSysPFMGSetup,(ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x));
2338:   return(0);
2339: }

2341: /*MC
2342:      PCSysPFMG - the hypre SysPFMG multigrid solver

2344:    Level: advanced

2346:    Options Database:
2347: + -pc_syspfmg_its <its> number of iterations of SysPFMG to use as preconditioner
2348: . -pc_syspfmg_num_pre_relax <steps> number of smoothing steps before coarse grid
2349: . -pc_syspfmg_num_post_relax <steps> number of smoothing steps after coarse grid
2350: . -pc_syspfmg_tol <tol> tolerance of SysPFMG
2351: . -pc_syspfmg_relax_type -relaxation type for the up and down cycles, one of Weighted-Jacobi,Red/Black-Gauss-Seidel

2353:    Notes:  This is for CELL-centered descretizations

2355:            This must be used with the MATHYPRESSTRUCT matrix type.
2356:            This is less general than in hypre, it supports only one part, and one block per process defined by a PETSc DMDA.
2357:            Also, only cell-centered variables.

2359: .seealso:  PCMG, MATHYPRESSTRUCT
2360: M*/

2362: PETSC_EXTERN PetscErrorCode PCCreate_SysPFMG(PC pc)
2363: {
2365:   PC_SysPFMG     *ex;

2368:   PetscNew(&ex); \
2369:   pc->data = ex;

2371:   ex->its            = 1;
2372:   ex->tol            = 1.e-8;
2373:   ex->relax_type     = 1;
2374:   ex->num_pre_relax  = 1;
2375:   ex->num_post_relax = 1;

2377:   pc->ops->setfromoptions  = PCSetFromOptions_SysPFMG;
2378:   pc->ops->view            = PCView_SysPFMG;
2379:   pc->ops->destroy         = PCDestroy_SysPFMG;
2380:   pc->ops->apply           = PCApply_SysPFMG;
2381:   pc->ops->applyrichardson = PCApplyRichardson_SysPFMG;
2382:   pc->ops->setup           = PCSetUp_SysPFMG;

2384:   MPI_Comm_dup(PetscObjectComm((PetscObject)pc),&(ex->hcomm));
2385:   PetscStackCallStandard(HYPRE_SStructSysPFMGCreate,(ex->hcomm,&ex->ss_solver));
2386:   return(0);
2387: }