Actual source code: hypre.c
petsc-3.11.4 2019-09-28
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: PetscBool nodal_relax;
77: PetscInt nodal_relax_levels;
79: PetscInt nodal_coarsening;
80: PetscInt nodal_coarsening_diag;
81: PetscInt vec_interp_variant;
82: PetscInt vec_interp_qmax;
83: PetscBool vec_interp_smooth;
84: PetscInt interp_refine;
86: HYPRE_IJVector *hmnull;
87: HYPRE_ParVector *phmnull; /* near null space passed to hypre */
88: PetscInt n_hmnull;
89: Vec hmnull_constant;
90: 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 */
92: /* options for AS (Auxiliary Space preconditioners) */
93: PetscInt as_print;
94: PetscInt as_max_iter;
95: PetscReal as_tol;
96: PetscInt as_relax_type;
97: PetscInt as_relax_times;
98: PetscReal as_relax_weight;
99: PetscReal as_omega;
100: PetscInt as_amg_alpha_opts[5]; /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for vector Poisson (AMS) or Curl problem (ADS) */
101: PetscReal as_amg_alpha_theta; /* AMG strength for vector Poisson (AMS) or Curl problem (ADS) */
102: PetscInt as_amg_beta_opts[5]; /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for scalar Poisson (AMS) or vector Poisson (ADS) */
103: PetscReal as_amg_beta_theta; /* AMG strength for scalar Poisson (AMS) or vector Poisson (ADS) */
104: PetscInt ams_cycle_type;
105: PetscInt ads_cycle_type;
107: /* additional data */
108: Mat G; /* MatHYPRE */
109: Mat C; /* MatHYPRE */
110: Mat alpha_Poisson; /* MatHYPRE */
111: Mat beta_Poisson; /* MatHYPRE */
113: /* extra information for AMS */
114: PetscInt dim; /* geometrical dimension */
115: HYPRE_IJVector coords[3];
116: HYPRE_IJVector constants[3];
117: Mat RT_PiFull, RT_Pi[3];
118: Mat ND_PiFull, ND_Pi[3];
119: PetscBool ams_beta_is_zero;
120: PetscBool ams_beta_is_zero_part;
121: PetscInt ams_proj_freq;
122: } PC_HYPRE;
124: PetscErrorCode PCHYPREGetSolver(PC pc,HYPRE_Solver *hsolver)
125: {
126: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
129: *hsolver = jac->hsolver;
130: return(0);
131: }
133: /* Resets (frees) Hypre's representation of the near null space */
134: static PetscErrorCode PCHYPREResetNearNullSpace_Private(PC pc)
135: {
136: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
137: PetscInt i;
139: PETSC_UNUSED PetscScalar *petscvecarray;
142: for (i=0; i<jac->n_hmnull; i++) {
143: VecHYPRE_ParVectorReplacePointer(jac->hmnull[i],jac->hmnull_hypre_data_array[i],petscvecarray);
144: PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->hmnull[i]));
145: }
146: PetscFree(jac->hmnull);
147: PetscFree(jac->hmnull_hypre_data_array);
148: PetscFree(jac->phmnull);
149: VecDestroy(&jac->hmnull_constant);
150: jac->n_hmnull = 0;
151: return(0);
152: }
154: static PetscErrorCode PCSetUp_HYPRE(PC pc)
155: {
156: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
157: Mat_HYPRE *hjac;
158: HYPRE_ParCSRMatrix hmat;
159: HYPRE_ParVector bv,xv;
160: PetscBool ishypre;
161: PetscErrorCode ierr;
164: if (!jac->hypre_type) {
165: PCHYPRESetType(pc,"boomeramg");
166: }
168: PetscObjectTypeCompare((PetscObject)pc->pmat,MATHYPRE,&ishypre);
169: if (!ishypre) {
170: MatDestroy(&jac->hpmat);
171: MatConvert(pc->pmat,MATHYPRE,MAT_INITIAL_MATRIX,&jac->hpmat);
172: } else {
173: PetscObjectReference((PetscObject)pc->pmat);
174: MatDestroy(&jac->hpmat);
175: jac->hpmat = pc->pmat;
176: }
177: hjac = (Mat_HYPRE*)(jac->hpmat->data);
179: /* special case for BoomerAMG */
180: if (jac->setup == HYPRE_BoomerAMGSetup) {
181: MatNullSpace mnull;
182: PetscBool has_const;
183: PetscInt bs,nvec,i;
184: const Vec *vecs;
185: PetscScalar *petscvecarray;
187: MatGetBlockSize(pc->pmat,&bs);
188: if (bs > 1) PetscStackCallStandard(HYPRE_BoomerAMGSetNumFunctions,(jac->hsolver,bs));
189: MatGetNearNullSpace(pc->mat, &mnull);
190: if (mnull) {
191: PCHYPREResetNearNullSpace_Private(pc);
192: MatNullSpaceGetVecs(mnull, &has_const, &nvec, &vecs);
193: PetscMalloc1(nvec+1,&jac->hmnull);
194: PetscMalloc1(nvec+1,&jac->hmnull_hypre_data_array);
195: PetscMalloc1(nvec+1,&jac->phmnull);
196: for (i=0; i<nvec; i++) {
197: VecHYPRE_IJVectorCreate(vecs[i],&jac->hmnull[i]);
198: VecGetArrayRead(vecs[i],(const PetscScalar **)&petscvecarray);
199: VecHYPRE_ParVectorReplacePointer(jac->hmnull[i],petscvecarray,jac->hmnull_hypre_data_array[i]);
200: VecRestoreArrayRead(vecs[i],(const PetscScalar **)&petscvecarray);
201: PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->hmnull[i],(void**)&jac->phmnull[i]));
202: }
203: if (has_const) {
204: MatCreateVecs(pc->pmat,&jac->hmnull_constant,NULL);
205: VecSet(jac->hmnull_constant,1);
206: VecNormalize(jac->hmnull_constant,NULL);
207: VecHYPRE_IJVectorCreate(jac->hmnull_constant,&jac->hmnull[nvec]);
208: VecGetArrayRead(jac->hmnull_constant,(const PetscScalar **)&petscvecarray);
209: VecHYPRE_ParVectorReplacePointer(jac->hmnull[nvec],petscvecarray,jac->hmnull_hypre_data_array[nvec]);
210: VecRestoreArrayRead(jac->hmnull_constant,(const PetscScalar **)&petscvecarray);
211: PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->hmnull[nvec],(void**)&jac->phmnull[nvec]));
212: nvec++;
213: }
214: PetscStackCallStandard(HYPRE_BoomerAMGSetInterpVectors,(jac->hsolver,nvec,jac->phmnull));
215: jac->n_hmnull = nvec;
216: }
217: }
219: /* special case for AMS */
220: if (jac->setup == HYPRE_AMSSetup) {
221: Mat_HYPRE *hm;
222: HYPRE_ParCSRMatrix parcsr;
223: if (!jac->coords[0] && !jac->constants[0] && !(jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1]))) {
224: 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");
225: }
226: if (jac->dim) {
227: PetscStackCallStandard(HYPRE_AMSSetDimension,(jac->hsolver,jac->dim));
228: }
229: if (jac->constants[0]) {
230: HYPRE_ParVector ozz,zoz,zzo = NULL;
231: PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->constants[0],(void**)(&ozz)));
232: PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->constants[1],(void**)(&zoz)));
233: if (jac->constants[2]) {
234: PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->constants[2],(void**)(&zzo)));
235: }
236: PetscStackCallStandard(HYPRE_AMSSetEdgeConstantVectors,(jac->hsolver,ozz,zoz,zzo));
237: }
238: if (jac->coords[0]) {
239: HYPRE_ParVector coords[3];
240: coords[0] = NULL;
241: coords[1] = NULL;
242: coords[2] = NULL;
243: if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[0],(void**)(&coords[0])));
244: if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[1],(void**)(&coords[1])));
245: if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[2],(void**)(&coords[2])));
246: PetscStackCallStandard(HYPRE_AMSSetCoordinateVectors,(jac->hsolver,coords[0],coords[1],coords[2]));
247: }
248: if (!jac->G) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE AMS preconditioner needs the discrete gradient operator via PCHYPRESetDiscreteGradient");
249: hm = (Mat_HYPRE*)(jac->G->data);
250: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&parcsr)));
251: PetscStackCallStandard(HYPRE_AMSSetDiscreteGradient,(jac->hsolver,parcsr));
252: if (jac->alpha_Poisson) {
253: hm = (Mat_HYPRE*)(jac->alpha_Poisson->data);
254: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&parcsr)));
255: PetscStackCallStandard(HYPRE_AMSSetAlphaPoissonMatrix,(jac->hsolver,parcsr));
256: }
257: if (jac->ams_beta_is_zero) {
258: PetscStackCallStandard(HYPRE_AMSSetBetaPoissonMatrix,(jac->hsolver,NULL));
259: } else if (jac->beta_Poisson) {
260: hm = (Mat_HYPRE*)(jac->beta_Poisson->data);
261: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&parcsr)));
262: PetscStackCallStandard(HYPRE_AMSSetBetaPoissonMatrix,(jac->hsolver,parcsr));
263: }
264: if (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1])) {
265: PetscInt i;
266: HYPRE_ParCSRMatrix nd_parcsrfull, nd_parcsr[3];
267: if (jac->ND_PiFull) {
268: hm = (Mat_HYPRE*)(jac->ND_PiFull->data);
269: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&nd_parcsrfull)));
270: } else {
271: nd_parcsrfull = NULL;
272: }
273: for (i=0;i<3;++i) {
274: if (jac->ND_Pi[i]) {
275: hm = (Mat_HYPRE*)(jac->ND_Pi[i]->data);
276: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&nd_parcsr[i])));
277: } else {
278: nd_parcsr[i] = NULL;
279: }
280: }
281: PetscStackCallStandard(HYPRE_AMSSetInterpolations,(jac->hsolver,nd_parcsrfull,nd_parcsr[0],nd_parcsr[1],nd_parcsr[2]));
282: }
283: }
284: /* special case for ADS */
285: if (jac->setup == HYPRE_ADSSetup) {
286: Mat_HYPRE *hm;
287: HYPRE_ParCSRMatrix parcsr;
288: 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])))) {
289: SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs either the coordinate vectors via PCSetCoordinates() or the interpolation matrices via PCHYPRESetInterpolations");
290: }
291: 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");
292: if (!jac->G) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs the discrete gradient operator via PCHYPRESetDiscreteGradient");
293: if (!jac->C) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs the discrete curl operator via PCHYPRESetDiscreteGradient");
294: if (jac->coords[0]) {
295: HYPRE_ParVector coords[3];
296: coords[0] = NULL;
297: coords[1] = NULL;
298: coords[2] = NULL;
299: if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[0],(void**)(&coords[0])));
300: if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[1],(void**)(&coords[1])));
301: if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[2],(void**)(&coords[2])));
302: PetscStackCallStandard(HYPRE_ADSSetCoordinateVectors,(jac->hsolver,coords[0],coords[1],coords[2]));
303: }
304: hm = (Mat_HYPRE*)(jac->G->data);
305: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&parcsr)));
306: PetscStackCallStandard(HYPRE_ADSSetDiscreteGradient,(jac->hsolver,parcsr));
307: hm = (Mat_HYPRE*)(jac->C->data);
308: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&parcsr)));
309: PetscStackCallStandard(HYPRE_ADSSetDiscreteCurl,(jac->hsolver,parcsr));
310: if ((jac->RT_PiFull || (jac->RT_Pi[0] && jac->RT_Pi[1])) && (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1]))) {
311: PetscInt i;
312: HYPRE_ParCSRMatrix rt_parcsrfull, rt_parcsr[3];
313: HYPRE_ParCSRMatrix nd_parcsrfull, nd_parcsr[3];
314: if (jac->RT_PiFull) {
315: hm = (Mat_HYPRE*)(jac->RT_PiFull->data);
316: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&rt_parcsrfull)));
317: } else {
318: rt_parcsrfull = NULL;
319: }
320: for (i=0;i<3;++i) {
321: if (jac->RT_Pi[i]) {
322: hm = (Mat_HYPRE*)(jac->RT_Pi[i]->data);
323: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&rt_parcsr[i])));
324: } else {
325: rt_parcsr[i] = NULL;
326: }
327: }
328: if (jac->ND_PiFull) {
329: hm = (Mat_HYPRE*)(jac->ND_PiFull->data);
330: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&nd_parcsrfull)));
331: } else {
332: nd_parcsrfull = NULL;
333: }
334: for (i=0;i<3;++i) {
335: if (jac->ND_Pi[i]) {
336: hm = (Mat_HYPRE*)(jac->ND_Pi[i]->data);
337: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hm->ij,(void**)(&nd_parcsr[i])));
338: } else {
339: nd_parcsr[i] = NULL;
340: }
341: }
342: 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]));
343: }
344: }
345: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hjac->ij,(void**)&hmat));
346: PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->b,(void**)&bv));
347: PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->x,(void**)&xv));
348: PetscStackCallStandard(jac->setup,(jac->hsolver,hmat,bv,xv));
349: return(0);
350: }
352: static PetscErrorCode PCApply_HYPRE(PC pc,Vec b,Vec x)
353: {
354: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
355: Mat_HYPRE *hjac = (Mat_HYPRE*)(jac->hpmat->data);
356: PetscErrorCode ierr;
357: HYPRE_ParCSRMatrix hmat;
358: PetscScalar *xv;
359: const PetscScalar *bv,*sbv;
360: HYPRE_ParVector jbv,jxv;
361: PetscScalar *sxv;
362: PetscInt hierr;
365: PetscCitationsRegister(hypreCitation,&cite);
366: if (!jac->applyrichardson) {VecSet(x,0.0);}
367: VecGetArrayRead(b,&bv);
368: VecGetArray(x,&xv);
369: VecHYPRE_ParVectorReplacePointer(hjac->b,(PetscScalar*)bv,sbv);
370: VecHYPRE_ParVectorReplacePointer(hjac->x,xv,sxv);
372: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hjac->ij,(void**)&hmat));
373: PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->b,(void**)&jbv));
374: PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->x,(void**)&jxv));
375: PetscStackCall("Hypre solve",h(*jac->solve)(jac->hsolver,hmat,jbv,jxv);
376: if (hierr && hierr != HYPRE_ERROR_CONV) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"Error in HYPRE solver, error code %d",hierr);
377: if (hierr) hypre__global_error = 0;);
379: if (jac->setup == HYPRE_AMSSetup && jac->ams_beta_is_zero_part) {
380: PetscStackCallStandard(HYPRE_AMSProjectOutGradients,(jac->hsolver,jxv));
381: }
382: VecHYPRE_ParVectorReplacePointer(hjac->b,(PetscScalar*)sbv,bv);
383: VecHYPRE_ParVectorReplacePointer(hjac->x,sxv,xv);
384: VecRestoreArray(x,&xv);
385: VecRestoreArrayRead(b,&bv);
386: return(0);
387: }
389: static PetscErrorCode PCReset_HYPRE(PC pc)
390: {
391: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
395: MatDestroy(&jac->hpmat);
396: MatDestroy(&jac->G);
397: MatDestroy(&jac->C);
398: MatDestroy(&jac->alpha_Poisson);
399: MatDestroy(&jac->beta_Poisson);
400: MatDestroy(&jac->RT_PiFull);
401: MatDestroy(&jac->RT_Pi[0]);
402: MatDestroy(&jac->RT_Pi[1]);
403: MatDestroy(&jac->RT_Pi[2]);
404: MatDestroy(&jac->ND_PiFull);
405: MatDestroy(&jac->ND_Pi[0]);
406: MatDestroy(&jac->ND_Pi[1]);
407: MatDestroy(&jac->ND_Pi[2]);
408: if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[0])); jac->coords[0] = NULL;
409: if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[1])); jac->coords[1] = NULL;
410: if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[2])); jac->coords[2] = NULL;
411: if (jac->constants[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[0])); jac->constants[0] = NULL;
412: if (jac->constants[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[1])); jac->constants[1] = NULL;
413: if (jac->constants[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[2])); jac->constants[2] = NULL;
414: PCHYPREResetNearNullSpace_Private(pc);
415: jac->ams_beta_is_zero = PETSC_FALSE;
416: jac->dim = 0;
417: return(0);
418: }
420: static PetscErrorCode PCDestroy_HYPRE(PC pc)
421: {
422: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
423: PetscErrorCode ierr;
426: PCReset_HYPRE(pc);
427: if (jac->destroy) PetscStackCallStandard(jac->destroy,(jac->hsolver));
428: PetscFree(jac->hypre_type);
429: if (jac->comm_hypre != MPI_COMM_NULL) { MPI_Comm_free(&(jac->comm_hypre));}
430: PetscFree(pc->data);
432: PetscObjectChangeTypeName((PetscObject)pc,0);
433: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetType_C",NULL);
434: PetscObjectComposeFunction((PetscObject)pc,"PCHYPREGetType_C",NULL);
435: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetCoordinates_C",NULL);
436: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteGradient_C",NULL);
437: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteCurl_C",NULL);
438: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetInterpolations_C",NULL);
439: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetConstantEdgeVectors_C",NULL);
440: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetPoissonMatrix_C",NULL);
441: return(0);
442: }
444: /* --------------------------------------------------------------------------------------------*/
445: static PetscErrorCode PCSetFromOptions_HYPRE_Pilut(PetscOptionItems *PetscOptionsObject,PC pc)
446: {
447: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
449: PetscBool flag;
452: PetscOptionsHead(PetscOptionsObject,"HYPRE Pilut Options");
453: PetscOptionsInt("-pc_hypre_pilut_maxiter","Number of iterations","None",jac->maxiter,&jac->maxiter,&flag);
454: if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetMaxIter,(jac->hsolver,jac->maxiter));
455: PetscOptionsReal("-pc_hypre_pilut_tol","Drop tolerance","None",jac->tol,&jac->tol,&flag);
456: if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetDropTolerance,(jac->hsolver,jac->tol));
457: PetscOptionsInt("-pc_hypre_pilut_factorrowsize","FactorRowSize","None",jac->factorrowsize,&jac->factorrowsize,&flag);
458: if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetFactorRowSize,(jac->hsolver,jac->factorrowsize));
459: PetscOptionsTail();
460: return(0);
461: }
463: static PetscErrorCode PCView_HYPRE_Pilut(PC pc,PetscViewer viewer)
464: {
465: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
467: PetscBool iascii;
470: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
471: if (iascii) {
472: PetscViewerASCIIPrintf(viewer," HYPRE Pilut preconditioning\n");
473: if (jac->maxiter != PETSC_DEFAULT) {
474: PetscViewerASCIIPrintf(viewer," maximum number of iterations %d\n",jac->maxiter);
475: } else {
476: PetscViewerASCIIPrintf(viewer," default maximum number of iterations \n");
477: }
478: if (jac->tol != PETSC_DEFAULT) {
479: PetscViewerASCIIPrintf(viewer," drop tolerance %g\n",(double)jac->tol);
480: } else {
481: PetscViewerASCIIPrintf(viewer," default drop tolerance \n");
482: }
483: if (jac->factorrowsize != PETSC_DEFAULT) {
484: PetscViewerASCIIPrintf(viewer," factor row size %d\n",jac->factorrowsize);
485: } else {
486: PetscViewerASCIIPrintf(viewer," default factor row size \n");
487: }
488: }
489: return(0);
490: }
492: /* --------------------------------------------------------------------------------------------*/
493: static PetscErrorCode PCSetFromOptions_HYPRE_Euclid(PetscOptionItems *PetscOptionsObject,PC pc)
494: {
495: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
497: PetscBool flag;
500: PetscOptionsHead(PetscOptionsObject,"HYPRE Euclid Options");
501: PetscOptionsInt("-pc_hypre_euclid_level","Factorization levels","None",jac->eu_level,&jac->eu_level,&flag);
502: if (flag) PetscStackCallStandard(HYPRE_EuclidSetLevel,(jac->hsolver,jac->eu_level));
503: PetscOptionsTail();
504: return(0);
505: }
507: static PetscErrorCode PCView_HYPRE_Euclid(PC pc,PetscViewer viewer)
508: {
509: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
511: PetscBool iascii;
514: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
515: if (iascii) {
516: PetscViewerASCIIPrintf(viewer," HYPRE Euclid preconditioning\n");
517: if (jac->eu_level != PETSC_DEFAULT) {
518: PetscViewerASCIIPrintf(viewer," factorization levels %d\n",jac->eu_level);
519: } else {
520: PetscViewerASCIIPrintf(viewer," default factorization levels \n");
521: }
522: }
523: return(0);
524: }
526: /* --------------------------------------------------------------------------------------------*/
528: static PetscErrorCode PCApplyTranspose_HYPRE_BoomerAMG(PC pc,Vec b,Vec x)
529: {
530: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
531: Mat_HYPRE *hjac = (Mat_HYPRE*)(jac->hpmat->data);
532: PetscErrorCode ierr;
533: HYPRE_ParCSRMatrix hmat;
534: PetscScalar *xv;
535: const PetscScalar *bv;
536: HYPRE_ParVector jbv,jxv;
537: PetscScalar *sbv,*sxv;
538: PetscInt hierr;
541: PetscCitationsRegister(hypreCitation,&cite);
542: VecSet(x,0.0);
543: VecGetArrayRead(b,&bv);
544: VecGetArray(x,&xv);
545: VecHYPRE_ParVectorReplacePointer(hjac->b,(PetscScalar*)bv,sbv);
546: VecHYPRE_ParVectorReplacePointer(hjac->x,xv,sxv);
548: PetscStackCallStandard(HYPRE_IJMatrixGetObject,(hjac->ij,(void**)&hmat));
549: PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->b,(void**)&jbv));
550: PetscStackCallStandard(HYPRE_IJVectorGetObject,(hjac->x,(void**)&jxv));
552: hHYPRE_BoomerAMGSolveT(jac->hsolver,hmat,jbv,jxv);
553: /* error code of 1 in BoomerAMG merely means convergence not achieved */
554: if (hierr && (hierr != 1)) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"Error in HYPRE solver, error code %d",hierr);
555: if (hierr) hypre__global_error = 0;
557: VecHYPRE_ParVectorReplacePointer(hjac->b,sbv,bv);
558: VecHYPRE_ParVectorReplacePointer(hjac->x,sxv,xv);
559: VecRestoreArray(x,&xv);
560: VecRestoreArrayRead(b,&bv);
561: return(0);
562: }
564: /* static array length */
565: #define ALEN(a) (sizeof(a)/sizeof((a)[0]))
567: static const char *HYPREBoomerAMGCycleType[] = {"","V","W"};
568: static const char *HYPREBoomerAMGCoarsenType[] = {"CLJP","Ruge-Stueben","","modifiedRuge-Stueben","","","Falgout", "", "PMIS", "", "HMIS"};
569: static const char *HYPREBoomerAMGMeasureType[] = {"local","global"};
570: /* The following corresponds to HYPRE_BoomerAMGSetRelaxType which has many missing numbers in the enum */
571: static const char *HYPREBoomerAMGSmoothType[] = {"Schwarz-smoothers","Pilut","ParaSails","Euclid"};
572: static const char *HYPREBoomerAMGRelaxType[] = {"Jacobi","sequential-Gauss-Seidel","seqboundary-Gauss-Seidel","SOR/Jacobi","backward-SOR/Jacobi",
573: "" /* [5] hybrid chaotic Gauss-Seidel (works only with OpenMP) */,"symmetric-SOR/Jacobi",
574: "" /* 7 */,"l1scaled-SOR/Jacobi","Gaussian-elimination",
575: "" /* 10 */, "" /* 11 */, "" /* 12 */, "l1-Gauss-Seidel" /* nonsymmetric */, "backward-l1-Gauss-Seidel" /* nonsymmetric */,
576: "CG" /* non-stationary */,"Chebyshev","FCF-Jacobi","l1scaled-Jacobi"};
577: static const char *HYPREBoomerAMGInterpType[] = {"classical", "", "", "direct", "multipass", "multipass-wts", "ext+i",
578: "ext+i-cc", "standard", "standard-wts", "block", "block-wtd", "FF", "FF1"};
579: static PetscErrorCode PCSetFromOptions_HYPRE_BoomerAMG(PetscOptionItems *PetscOptionsObject,PC pc)
580: {
581: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
583: PetscInt bs,n,indx,level;
584: PetscBool flg, tmp_truth;
585: double tmpdbl, twodbl[2];
588: PetscOptionsHead(PetscOptionsObject,"HYPRE BoomerAMG Options");
589: PetscOptionsEList("-pc_hypre_boomeramg_cycle_type","Cycle type","None",HYPREBoomerAMGCycleType+1,2,HYPREBoomerAMGCycleType[jac->cycletype],&indx,&flg);
590: if (flg) {
591: jac->cycletype = indx+1;
592: PetscStackCallStandard(HYPRE_BoomerAMGSetCycleType,(jac->hsolver,jac->cycletype));
593: }
594: PetscOptionsInt("-pc_hypre_boomeramg_max_levels","Number of levels (of grids) allowed","None",jac->maxlevels,&jac->maxlevels,&flg);
595: if (flg) {
596: if (jac->maxlevels < 2) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of levels %d must be at least two",jac->maxlevels);
597: PetscStackCallStandard(HYPRE_BoomerAMGSetMaxLevels,(jac->hsolver,jac->maxlevels));
598: }
599: PetscOptionsInt("-pc_hypre_boomeramg_max_iter","Maximum iterations used PER hypre call","None",jac->maxiter,&jac->maxiter,&flg);
600: if (flg) {
601: if (jac->maxiter < 1) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of iterations %d must be at least one",jac->maxiter);
602: PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,jac->maxiter));
603: }
604: PetscOptionsScalar("-pc_hypre_boomeramg_tol","Convergence tolerance PER hypre call (0.0 = use a fixed number of iterations)","None",jac->tol,&jac->tol,&flg);
605: if (flg) {
606: 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);
607: PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,jac->tol));
608: }
609: bs = 1;
610: if (pc->pmat) {
611: MatGetBlockSize(pc->pmat,&bs);
612: }
613: PetscOptionsInt("-pc_hypre_boomeramg_numfunctions","Number of functions","HYPRE_BoomerAMGSetNumFunctions",bs,&bs,&flg);
614: if (flg) {
615: PetscStackCallStandard(HYPRE_BoomerAMGSetNumFunctions,(jac->hsolver,bs));
616: }
618: PetscOptionsScalar("-pc_hypre_boomeramg_truncfactor","Truncation factor for interpolation (0=no truncation)","None",jac->truncfactor,&jac->truncfactor,&flg);
619: if (flg) {
620: 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);
621: PetscStackCallStandard(HYPRE_BoomerAMGSetTruncFactor,(jac->hsolver,jac->truncfactor));
622: }
624: PetscOptionsInt("-pc_hypre_boomeramg_P_max","Max elements per row for interpolation operator (0=unlimited)","None",jac->pmax,&jac->pmax,&flg);
625: if (flg) {
626: 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);
627: PetscStackCallStandard(HYPRE_BoomerAMGSetPMaxElmts,(jac->hsolver,jac->pmax));
628: }
630: PetscOptionsInt("-pc_hypre_boomeramg_agg_nl","Number of levels of aggressive coarsening","None",jac->agg_nl,&jac->agg_nl,&flg);
631: if (flg) {
632: 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);
634: PetscStackCallStandard(HYPRE_BoomerAMGSetAggNumLevels,(jac->hsolver,jac->agg_nl));
635: }
637: PetscOptionsInt("-pc_hypre_boomeramg_agg_num_paths","Number of paths for aggressive coarsening","None",jac->agg_num_paths,&jac->agg_num_paths,&flg);
638: if (flg) {
639: 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);
641: PetscStackCallStandard(HYPRE_BoomerAMGSetNumPaths,(jac->hsolver,jac->agg_num_paths));
642: }
645: PetscOptionsScalar("-pc_hypre_boomeramg_strong_threshold","Threshold for being strongly connected","None",jac->strongthreshold,&jac->strongthreshold,&flg);
646: if (flg) {
647: 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);
648: PetscStackCallStandard(HYPRE_BoomerAMGSetStrongThreshold,(jac->hsolver,jac->strongthreshold));
649: }
650: PetscOptionsScalar("-pc_hypre_boomeramg_max_row_sum","Maximum row sum","None",jac->maxrowsum,&jac->maxrowsum,&flg);
651: if (flg) {
652: 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);
653: 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);
654: PetscStackCallStandard(HYPRE_BoomerAMGSetMaxRowSum,(jac->hsolver,jac->maxrowsum));
655: }
657: /* Grid sweeps */
658: PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_all","Number of sweeps for the up and down grid levels","None",jac->gridsweeps[0],&indx,&flg);
659: if (flg) {
660: PetscStackCallStandard(HYPRE_BoomerAMGSetNumSweeps,(jac->hsolver,indx));
661: /* modify the jac structure so we can view the updated options with PC_View */
662: jac->gridsweeps[0] = indx;
663: jac->gridsweeps[1] = indx;
664: /*defaults coarse to 1 */
665: jac->gridsweeps[2] = 1;
666: }
667: PetscOptionsInt("-pc_hypre_boomeramg_nodal_coarsen","Use a nodal based coarsening 1-6","HYPRE_BoomerAMGSetNodal",jac->nodal_coarsening,&jac->nodal_coarsening,&flg);
668: if (flg) {
669: PetscStackCallStandard(HYPRE_BoomerAMGSetNodal,(jac->hsolver,jac->nodal_coarsening));
670: }
671: PetscOptionsInt("-pc_hypre_boomeramg_nodal_coarsen_diag","Diagonal in strength matrix for nodal based coarsening 0-2","HYPRE_BoomerAMGSetNodalDiag",jac->nodal_coarsening_diag,&jac->nodal_coarsening_diag,&flg);
672: if (flg) {
673: PetscStackCallStandard(HYPRE_BoomerAMGSetNodalDiag,(jac->hsolver,jac->nodal_coarsening_diag));
674: }
675: PetscOptionsInt("-pc_hypre_boomeramg_vec_interp_variant","Variant of algorithm 1-3","HYPRE_BoomerAMGSetInterpVecVariant",jac->vec_interp_variant, &jac->vec_interp_variant,&flg);
676: if (flg) {
677: PetscStackCallStandard(HYPRE_BoomerAMGSetInterpVecVariant,(jac->hsolver,jac->vec_interp_variant));
678: }
679: PetscOptionsInt("-pc_hypre_boomeramg_vec_interp_qmax","Max elements per row for each Q","HYPRE_BoomerAMGSetInterpVecQMax",jac->vec_interp_qmax, &jac->vec_interp_qmax,&flg);
680: if (flg) {
681: PetscStackCallStandard(HYPRE_BoomerAMGSetInterpVecQMax,(jac->hsolver,jac->vec_interp_qmax));
682: }
683: PetscOptionsBool("-pc_hypre_boomeramg_vec_interp_smooth","Whether to smooth the interpolation vectors","HYPRE_BoomerAMGSetSmoothInterpVectors",jac->vec_interp_smooth, &jac->vec_interp_smooth,&flg);
684: if (flg) {
685: PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothInterpVectors,(jac->hsolver,jac->vec_interp_smooth));
686: }
687: PetscOptionsInt("-pc_hypre_boomeramg_interp_refine","Preprocess the interpolation matrix through iterative weight refinement","HYPRE_BoomerAMGSetInterpRefine",jac->interp_refine, &jac->interp_refine,&flg);
688: if (flg) {
689: PetscStackCallStandard(HYPRE_BoomerAMGSetInterpRefine,(jac->hsolver,jac->interp_refine));
690: }
691: PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_down","Number of sweeps for the down cycles","None",jac->gridsweeps[0], &indx,&flg);
692: if (flg) {
693: PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,(jac->hsolver,indx, 1));
694: jac->gridsweeps[0] = indx;
695: }
696: PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_up","Number of sweeps for the up cycles","None",jac->gridsweeps[1],&indx,&flg);
697: if (flg) {
698: PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,(jac->hsolver,indx, 2));
699: jac->gridsweeps[1] = indx;
700: }
701: PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_coarse","Number of sweeps for the coarse level","None",jac->gridsweeps[2],&indx,&flg);
702: if (flg) {
703: PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,(jac->hsolver,indx, 3));
704: jac->gridsweeps[2] = indx;
705: }
707: /* Smooth type */
708: PetscOptionsEList("-pc_hypre_boomeramg_smooth_type","Enable more complex smoothers","None",HYPREBoomerAMGSmoothType,ALEN(HYPREBoomerAMGSmoothType),HYPREBoomerAMGSmoothType[0],&indx,&flg);
709: if (flg) {
710: jac->smoothtype = indx;
711: PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothType,(jac->hsolver,indx+6));
712: jac->smoothnumlevels = 25;
713: PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothNumLevels,(jac->hsolver,25));
714: }
716: /* Number of smoothing levels */
717: PetscOptionsInt("-pc_hypre_boomeramg_smooth_num_levels","Number of levels on which more complex smoothers are used","None",25,&indx,&flg);
718: if (flg && (jac->smoothtype != -1)) {
719: jac->smoothnumlevels = indx;
720: PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothNumLevels,(jac->hsolver,indx));
721: }
723: /* Number of levels for ILU(k) for Euclid */
724: PetscOptionsInt("-pc_hypre_boomeramg_eu_level","Number of levels for ILU(k) in Euclid smoother","None",0,&indx,&flg);
725: if (flg && (jac->smoothtype == 3)) {
726: jac->eu_level = indx;
727: PetscStackCallStandard(HYPRE_BoomerAMGSetEuLevel,(jac->hsolver,indx));
728: }
730: /* Filter for ILU(k) for Euclid */
731: double droptolerance;
732: PetscOptionsScalar("-pc_hypre_boomeramg_eu_droptolerance","Drop tolerance for ILU(k) in Euclid smoother","None",0,&droptolerance,&flg);
733: if (flg && (jac->smoothtype == 3)) {
734: jac->eu_droptolerance = droptolerance;
735: PetscStackCallStandard(HYPRE_BoomerAMGSetEuLevel,(jac->hsolver,droptolerance));
736: }
738: /* Use Block Jacobi ILUT for Euclid */
739: PetscOptionsBool("-pc_hypre_boomeramg_eu_bj", "Use Block Jacobi for ILU in Euclid smoother?", "None", PETSC_FALSE, &tmp_truth, &flg);
740: if (flg && (jac->smoothtype == 3)) {
741: jac->eu_bj = tmp_truth;
742: PetscStackCallStandard(HYPRE_BoomerAMGSetEuBJ,(jac->hsolver,jac->eu_bj));
743: }
745: /* Relax type */
746: PetscOptionsEList("-pc_hypre_boomeramg_relax_type_all","Relax type for the up and down cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg);
747: if (flg) {
748: jac->relaxtype[0] = jac->relaxtype[1] = indx;
749: PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxType,(jac->hsolver, indx));
750: /* by default, coarse type set to 9 */
751: jac->relaxtype[2] = 9;
752: PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,(jac->hsolver, 9, 3));
753: }
754: PetscOptionsEList("-pc_hypre_boomeramg_relax_type_down","Relax type for the down cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg);
755: if (flg) {
756: jac->relaxtype[0] = indx;
757: PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,(jac->hsolver, indx, 1));
758: }
759: PetscOptionsEList("-pc_hypre_boomeramg_relax_type_up","Relax type for the up cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg);
760: if (flg) {
761: jac->relaxtype[1] = indx;
762: PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,(jac->hsolver, indx, 2));
763: }
764: PetscOptionsEList("-pc_hypre_boomeramg_relax_type_coarse","Relax type on coarse grid","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[9],&indx,&flg);
765: if (flg) {
766: jac->relaxtype[2] = indx;
767: PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,(jac->hsolver, indx, 3));
768: }
770: /* Relaxation Weight */
771: 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);
772: if (flg) {
773: PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxWt,(jac->hsolver,tmpdbl));
774: jac->relaxweight = tmpdbl;
775: }
777: n = 2;
778: twodbl[0] = twodbl[1] = 1.0;
779: PetscOptionsRealArray("-pc_hypre_boomeramg_relax_weight_level","Set the relaxation weight for a particular level (weight,level)","None",twodbl, &n, &flg);
780: if (flg) {
781: if (n == 2) {
782: indx = (int)PetscAbsReal(twodbl[1]);
783: PetscStackCallStandard(HYPRE_BoomerAMGSetLevelRelaxWt,(jac->hsolver,twodbl[0],indx));
784: } 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);
785: }
787: /* Outer relaxation Weight */
788: 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);
789: if (flg) {
790: PetscStackCallStandard(HYPRE_BoomerAMGSetOuterWt,(jac->hsolver, tmpdbl));
791: jac->outerrelaxweight = tmpdbl;
792: }
794: n = 2;
795: twodbl[0] = twodbl[1] = 1.0;
796: PetscOptionsRealArray("-pc_hypre_boomeramg_outer_relax_weight_level","Set the outer relaxation weight for a particular level (weight,level)","None",twodbl, &n, &flg);
797: if (flg) {
798: if (n == 2) {
799: indx = (int)PetscAbsReal(twodbl[1]);
800: PetscStackCallStandard(HYPRE_BoomerAMGSetLevelOuterWt,(jac->hsolver, twodbl[0], indx));
801: } 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);
802: }
804: /* the Relax Order */
805: PetscOptionsBool("-pc_hypre_boomeramg_no_CF", "Do not use CF-relaxation", "None", PETSC_FALSE, &tmp_truth, &flg);
807: if (flg && tmp_truth) {
808: jac->relaxorder = 0;
809: PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxOrder,(jac->hsolver, jac->relaxorder));
810: }
811: PetscOptionsEList("-pc_hypre_boomeramg_measure_type","Measure type","None",HYPREBoomerAMGMeasureType,ALEN(HYPREBoomerAMGMeasureType),HYPREBoomerAMGMeasureType[0],&indx,&flg);
812: if (flg) {
813: jac->measuretype = indx;
814: PetscStackCallStandard(HYPRE_BoomerAMGSetMeasureType,(jac->hsolver,jac->measuretype));
815: }
816: /* update list length 3/07 */
817: PetscOptionsEList("-pc_hypre_boomeramg_coarsen_type","Coarsen type","None",HYPREBoomerAMGCoarsenType,ALEN(HYPREBoomerAMGCoarsenType),HYPREBoomerAMGCoarsenType[6],&indx,&flg);
818: if (flg) {
819: jac->coarsentype = indx;
820: PetscStackCallStandard(HYPRE_BoomerAMGSetCoarsenType,(jac->hsolver,jac->coarsentype));
821: }
823: /* new 3/07 */
824: PetscOptionsEList("-pc_hypre_boomeramg_interp_type","Interpolation type","None",HYPREBoomerAMGInterpType,ALEN(HYPREBoomerAMGInterpType),HYPREBoomerAMGInterpType[0],&indx,&flg);
825: if (flg) {
826: jac->interptype = indx;
827: PetscStackCallStandard(HYPRE_BoomerAMGSetInterpType,(jac->hsolver,jac->interptype));
828: }
830: PetscOptionsName("-pc_hypre_boomeramg_print_statistics","Print statistics","None",&flg);
831: if (flg) {
832: level = 3;
833: PetscOptionsInt("-pc_hypre_boomeramg_print_statistics","Print statistics","None",level,&level,NULL);
835: jac->printstatistics = PETSC_TRUE;
836: PetscStackCallStandard(HYPRE_BoomerAMGSetPrintLevel,(jac->hsolver,level));
837: }
839: PetscOptionsName("-pc_hypre_boomeramg_print_debug","Print debug information","None",&flg);
840: if (flg) {
841: level = 3;
842: PetscOptionsInt("-pc_hypre_boomeramg_print_debug","Print debug information","None",level,&level,NULL);
844: jac->printstatistics = PETSC_TRUE;
845: PetscStackCallStandard(HYPRE_BoomerAMGSetDebugFlag,(jac->hsolver,level));
846: }
848: PetscOptionsBool("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None", PETSC_FALSE, &tmp_truth, &flg);
849: if (flg && tmp_truth) {
850: PetscInt tmp_int;
851: PetscOptionsInt("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None",jac->nodal_relax_levels,&tmp_int,&flg);
852: if (flg) jac->nodal_relax_levels = tmp_int;
853: PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothType,(jac->hsolver,6));
854: PetscStackCallStandard(HYPRE_BoomerAMGSetDomainType,(jac->hsolver,1));
855: PetscStackCallStandard(HYPRE_BoomerAMGSetOverlap,(jac->hsolver,0));
856: PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothNumLevels,(jac->hsolver,jac->nodal_relax_levels));
857: }
859: PetscOptionsTail();
860: return(0);
861: }
863: 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)
864: {
865: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
867: PetscInt oits;
870: PetscCitationsRegister(hypreCitation,&cite);
871: PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,its*jac->maxiter));
872: PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,rtol));
873: jac->applyrichardson = PETSC_TRUE;
874: PCApply_HYPRE(pc,b,y);
875: jac->applyrichardson = PETSC_FALSE;
876: PetscStackCallStandard(HYPRE_BoomerAMGGetNumIterations,(jac->hsolver,(HYPRE_Int *)&oits));
877: *outits = oits;
878: if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
879: else *reason = PCRICHARDSON_CONVERGED_RTOL;
880: PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,jac->tol));
881: PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,jac->maxiter));
882: return(0);
883: }
886: static PetscErrorCode PCView_HYPRE_BoomerAMG(PC pc,PetscViewer viewer)
887: {
888: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
890: PetscBool iascii;
893: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
894: if (iascii) {
895: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG preconditioning\n");
896: PetscViewerASCIIPrintf(viewer," Cycle type %s\n",HYPREBoomerAMGCycleType[jac->cycletype]);
897: PetscViewerASCIIPrintf(viewer," Maximum number of levels %D\n",jac->maxlevels);
898: PetscViewerASCIIPrintf(viewer," Maximum number of iterations PER hypre call %D\n",jac->maxiter);
899: PetscViewerASCIIPrintf(viewer," Convergence tolerance PER hypre call %g\n",(double)jac->tol);
900: PetscViewerASCIIPrintf(viewer," Threshold for strong coupling %g\n",(double)jac->strongthreshold);
901: PetscViewerASCIIPrintf(viewer," Interpolation truncation factor %g\n",(double)jac->truncfactor);
902: PetscViewerASCIIPrintf(viewer," Interpolation: max elements per row %D\n",jac->pmax);
903: if (jac->interp_refine) {
904: PetscViewerASCIIPrintf(viewer," Interpolation: number of steps of weighted refinement %D\n",jac->interp_refine);
905: }
906: PetscViewerASCIIPrintf(viewer," Number of levels of aggressive coarsening %D\n",jac->agg_nl);
907: PetscViewerASCIIPrintf(viewer," Number of paths for aggressive coarsening %D\n",jac->agg_num_paths);
909: PetscViewerASCIIPrintf(viewer," Maximum row sums %g\n",(double)jac->maxrowsum);
911: PetscViewerASCIIPrintf(viewer," Sweeps down %D\n",jac->gridsweeps[0]);
912: PetscViewerASCIIPrintf(viewer," Sweeps up %D\n",jac->gridsweeps[1]);
913: PetscViewerASCIIPrintf(viewer," Sweeps on coarse %D\n",jac->gridsweeps[2]);
915: PetscViewerASCIIPrintf(viewer," Relax down %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[0]]);
916: PetscViewerASCIIPrintf(viewer," Relax up %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[1]]);
917: PetscViewerASCIIPrintf(viewer," Relax on coarse %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[2]]);
919: PetscViewerASCIIPrintf(viewer," Relax weight (all) %g\n",(double)jac->relaxweight);
920: PetscViewerASCIIPrintf(viewer," Outer relax weight (all) %g\n",(double)jac->outerrelaxweight);
922: if (jac->relaxorder) {
923: PetscViewerASCIIPrintf(viewer," Using CF-relaxation\n");
924: } else {
925: PetscViewerASCIIPrintf(viewer," Not using CF-relaxation\n");
926: }
927: if (jac->smoothtype!=-1) {
928: PetscViewerASCIIPrintf(viewer," Smooth type %s\n",HYPREBoomerAMGSmoothType[jac->smoothtype]);
929: PetscViewerASCIIPrintf(viewer," Smooth num levels %D\n",jac->smoothnumlevels);
930: } else {
931: PetscViewerASCIIPrintf(viewer," Not using more complex smoothers.\n");
932: }
933: if (jac->smoothtype==3) {
934: PetscViewerASCIIPrintf(viewer," Euclid ILU(k) levels %D\n",jac->eu_level);
935: PetscViewerASCIIPrintf(viewer," Euclid ILU(k) drop tolerance %g\n",(double)jac->eu_droptolerance);
936: PetscViewerASCIIPrintf(viewer," Euclid ILU use Block-Jacobi? %D\n",jac->eu_bj);
937: }
938: PetscViewerASCIIPrintf(viewer," Measure type %s\n",HYPREBoomerAMGMeasureType[jac->measuretype]);
939: PetscViewerASCIIPrintf(viewer," Coarsen type %s\n",HYPREBoomerAMGCoarsenType[jac->coarsentype]);
940: PetscViewerASCIIPrintf(viewer," Interpolation type %s\n",HYPREBoomerAMGInterpType[jac->interptype]);
941: if (jac->nodal_coarsening) {
942: PetscViewerASCIIPrintf(viewer," Using nodal coarsening (with HYPRE_BOOMERAMGSetNodal() %D\n",jac->nodal_coarsening);
943: }
944: if (jac->vec_interp_variant) {
945: PetscViewerASCIIPrintf(viewer," HYPRE_BoomerAMGSetInterpVecVariant() %D\n",jac->vec_interp_variant);
946: PetscViewerASCIIPrintf(viewer," HYPRE_BoomerAMGSetInterpVecQMax() %D\n",jac->vec_interp_qmax);
947: PetscViewerASCIIPrintf(viewer," HYPRE_BoomerAMGSetSmoothInterpVectors() %d\n",jac->vec_interp_smooth);
948: }
949: if (jac->nodal_relax) {
950: PetscViewerASCIIPrintf(viewer," Using nodal relaxation via Schwarz smoothing on levels %D\n",jac->nodal_relax_levels);
951: }
952: }
953: return(0);
954: }
956: /* --------------------------------------------------------------------------------------------*/
957: static PetscErrorCode PCSetFromOptions_HYPRE_ParaSails(PetscOptionItems *PetscOptionsObject,PC pc)
958: {
959: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
961: PetscInt indx;
962: PetscBool flag;
963: const char *symtlist[] = {"nonsymmetric","SPD","nonsymmetric,SPD"};
966: PetscOptionsHead(PetscOptionsObject,"HYPRE ParaSails Options");
967: PetscOptionsInt("-pc_hypre_parasails_nlevels","Number of number of levels","None",jac->nlevels,&jac->nlevels,0);
968: PetscOptionsReal("-pc_hypre_parasails_thresh","Threshold","None",jac->threshhold,&jac->threshhold,&flag);
969: if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetParams,(jac->hsolver,jac->threshhold,jac->nlevels));
971: PetscOptionsReal("-pc_hypre_parasails_filter","filter","None",jac->filter,&jac->filter,&flag);
972: if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetFilter,(jac->hsolver,jac->filter));
974: PetscOptionsReal("-pc_hypre_parasails_loadbal","Load balance","None",jac->loadbal,&jac->loadbal,&flag);
975: if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetLoadbal,(jac->hsolver,jac->loadbal));
977: PetscOptionsBool("-pc_hypre_parasails_logging","Print info to screen","None",(PetscBool)jac->logging,(PetscBool*)&jac->logging,&flag);
978: if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetLogging,(jac->hsolver,jac->logging));
980: PetscOptionsBool("-pc_hypre_parasails_reuse","Reuse nonzero pattern in preconditioner","None",(PetscBool)jac->ruse,(PetscBool*)&jac->ruse,&flag);
981: if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetReuse,(jac->hsolver,jac->ruse));
983: PetscOptionsEList("-pc_hypre_parasails_sym","Symmetry of matrix and preconditioner","None",symtlist,ALEN(symtlist),symtlist[0],&indx,&flag);
984: if (flag) {
985: jac->symt = indx;
986: PetscStackCallStandard(HYPRE_ParaSailsSetSym,(jac->hsolver,jac->symt));
987: }
989: PetscOptionsTail();
990: return(0);
991: }
993: static PetscErrorCode PCView_HYPRE_ParaSails(PC pc,PetscViewer viewer)
994: {
995: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
997: PetscBool iascii;
998: const char *symt = 0;;
1001: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
1002: if (iascii) {
1003: PetscViewerASCIIPrintf(viewer," HYPRE ParaSails preconditioning\n");
1004: PetscViewerASCIIPrintf(viewer," nlevels %d\n",jac->nlevels);
1005: PetscViewerASCIIPrintf(viewer," threshold %g\n",(double)jac->threshhold);
1006: PetscViewerASCIIPrintf(viewer," filter %g\n",(double)jac->filter);
1007: PetscViewerASCIIPrintf(viewer," load balance %g\n",(double)jac->loadbal);
1008: PetscViewerASCIIPrintf(viewer," reuse nonzero structure %s\n",PetscBools[jac->ruse]);
1009: PetscViewerASCIIPrintf(viewer," print info to screen %s\n",PetscBools[jac->logging]);
1010: if (!jac->symt) symt = "nonsymmetric matrix and preconditioner";
1011: else if (jac->symt == 1) symt = "SPD matrix and preconditioner";
1012: else if (jac->symt == 2) symt = "nonsymmetric matrix but SPD preconditioner";
1013: else SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG,"Unknown HYPRE ParaSails symmetric option %d",jac->symt);
1014: PetscViewerASCIIPrintf(viewer," %s\n",symt);
1015: }
1016: return(0);
1017: }
1018: /* --------------------------------------------------------------------------------------------*/
1019: static PetscErrorCode PCSetFromOptions_HYPRE_AMS(PetscOptionItems *PetscOptionsObject,PC pc)
1020: {
1021: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1023: PetscInt n;
1024: PetscBool flag,flag2,flag3,flag4;
1027: PetscOptionsHead(PetscOptionsObject,"HYPRE AMS Options");
1028: PetscOptionsInt("-pc_hypre_ams_print_level","Debugging output level for AMS","None",jac->as_print,&jac->as_print,&flag);
1029: if (flag) PetscStackCallStandard(HYPRE_AMSSetPrintLevel,(jac->hsolver,jac->as_print));
1030: PetscOptionsInt("-pc_hypre_ams_max_iter","Maximum number of AMS multigrid iterations within PCApply","None",jac->as_max_iter,&jac->as_max_iter,&flag);
1031: if (flag) PetscStackCallStandard(HYPRE_AMSSetMaxIter,(jac->hsolver,jac->as_max_iter));
1032: PetscOptionsInt("-pc_hypre_ams_cycle_type","Cycle type for AMS multigrid","None",jac->ams_cycle_type,&jac->ams_cycle_type,&flag);
1033: if (flag) PetscStackCallStandard(HYPRE_AMSSetCycleType,(jac->hsolver,jac->ams_cycle_type));
1034: PetscOptionsReal("-pc_hypre_ams_tol","Error tolerance for AMS multigrid","None",jac->as_tol,&jac->as_tol,&flag);
1035: if (flag) PetscStackCallStandard(HYPRE_AMSSetTol,(jac->hsolver,jac->as_tol));
1036: PetscOptionsInt("-pc_hypre_ams_relax_type","Relaxation type for AMS smoother","None",jac->as_relax_type,&jac->as_relax_type,&flag);
1037: PetscOptionsInt("-pc_hypre_ams_relax_times","Number of relaxation steps for AMS smoother","None",jac->as_relax_times,&jac->as_relax_times,&flag2);
1038: PetscOptionsReal("-pc_hypre_ams_relax_weight","Relaxation weight for AMS smoother","None",jac->as_relax_weight,&jac->as_relax_weight,&flag3);
1039: PetscOptionsReal("-pc_hypre_ams_omega","SSOR coefficient for AMS smoother","None",jac->as_omega,&jac->as_omega,&flag4);
1040: if (flag || flag2 || flag3 || flag4) {
1041: PetscStackCallStandard(HYPRE_AMSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type,
1042: jac->as_relax_times,
1043: jac->as_relax_weight,
1044: jac->as_omega));
1045: }
1046: 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);
1047: n = 5;
1048: PetscOptionsIntArray("-pc_hypre_ams_amg_alpha_options","AMG options for vector Poisson","None",jac->as_amg_alpha_opts,&n,&flag2);
1049: if (flag || flag2) {
1050: PetscStackCallStandard(HYPRE_AMSSetAlphaAMGOptions,(jac->hsolver,jac->as_amg_alpha_opts[0], /* AMG coarsen type */
1051: jac->as_amg_alpha_opts[1], /* AMG agg_levels */
1052: jac->as_amg_alpha_opts[2], /* AMG relax_type */
1053: jac->as_amg_alpha_theta,
1054: jac->as_amg_alpha_opts[3], /* AMG interp_type */
1055: jac->as_amg_alpha_opts[4])); /* AMG Pmax */
1056: }
1057: 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);
1058: n = 5;
1059: PetscOptionsIntArray("-pc_hypre_ams_amg_beta_options","AMG options for scalar Poisson solver","None",jac->as_amg_beta_opts,&n,&flag2);
1060: if (flag || flag2) {
1061: PetscStackCallStandard(HYPRE_AMSSetBetaAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */
1062: jac->as_amg_beta_opts[1], /* AMG agg_levels */
1063: jac->as_amg_beta_opts[2], /* AMG relax_type */
1064: jac->as_amg_beta_theta,
1065: jac->as_amg_beta_opts[3], /* AMG interp_type */
1066: jac->as_amg_beta_opts[4])); /* AMG Pmax */
1067: }
1068: 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);
1069: if (flag) { /* override HYPRE's default only if the options is used */
1070: PetscStackCallStandard(HYPRE_AMSSetProjectionFrequency,(jac->hsolver,jac->ams_proj_freq));
1071: }
1072: PetscOptionsTail();
1073: return(0);
1074: }
1076: static PetscErrorCode PCView_HYPRE_AMS(PC pc,PetscViewer viewer)
1077: {
1078: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1080: PetscBool iascii;
1083: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
1084: if (iascii) {
1085: PetscViewerASCIIPrintf(viewer," HYPRE AMS preconditioning\n");
1086: PetscViewerASCIIPrintf(viewer," subspace iterations per application %d\n",jac->as_max_iter);
1087: PetscViewerASCIIPrintf(viewer," subspace cycle type %d\n",jac->ams_cycle_type);
1088: PetscViewerASCIIPrintf(viewer," subspace iteration tolerance %g\n",jac->as_tol);
1089: PetscViewerASCIIPrintf(viewer," smoother type %d\n",jac->as_relax_type);
1090: PetscViewerASCIIPrintf(viewer," number of smoothing steps %d\n",jac->as_relax_times);
1091: PetscViewerASCIIPrintf(viewer," smoother weight %g\n",jac->as_relax_weight);
1092: PetscViewerASCIIPrintf(viewer," smoother omega %g\n",jac->as_omega);
1093: if (jac->alpha_Poisson) {
1094: PetscViewerASCIIPrintf(viewer," vector Poisson solver (passed in by user)\n");
1095: } else {
1096: PetscViewerASCIIPrintf(viewer," vector Poisson solver (computed) \n");
1097: }
1098: PetscViewerASCIIPrintf(viewer," boomerAMG coarsening type %d\n",jac->as_amg_alpha_opts[0]);
1099: PetscViewerASCIIPrintf(viewer," boomerAMG levels of aggressive coarsening %d\n",jac->as_amg_alpha_opts[1]);
1100: PetscViewerASCIIPrintf(viewer," boomerAMG relaxation type %d\n",jac->as_amg_alpha_opts[2]);
1101: PetscViewerASCIIPrintf(viewer," boomerAMG interpolation type %d\n",jac->as_amg_alpha_opts[3]);
1102: PetscViewerASCIIPrintf(viewer," boomerAMG max nonzero elements in interpolation rows %d\n",jac->as_amg_alpha_opts[4]);
1103: PetscViewerASCIIPrintf(viewer," boomerAMG strength threshold %g\n",jac->as_amg_alpha_theta);
1104: if (!jac->ams_beta_is_zero) {
1105: if (jac->beta_Poisson) {
1106: PetscViewerASCIIPrintf(viewer," scalar Poisson solver (passed in by user)\n");
1107: } else {
1108: PetscViewerASCIIPrintf(viewer," scalar Poisson solver (computed) \n");
1109: }
1110: PetscViewerASCIIPrintf(viewer," boomerAMG coarsening type %d\n",jac->as_amg_beta_opts[0]);
1111: PetscViewerASCIIPrintf(viewer," boomerAMG levels of aggressive coarsening %d\n",jac->as_amg_beta_opts[1]);
1112: PetscViewerASCIIPrintf(viewer," boomerAMG relaxation type %d\n",jac->as_amg_beta_opts[2]);
1113: PetscViewerASCIIPrintf(viewer," boomerAMG interpolation type %d\n",jac->as_amg_beta_opts[3]);
1114: PetscViewerASCIIPrintf(viewer," boomerAMG max nonzero elements in interpolation rows %d\n",jac->as_amg_beta_opts[4]);
1115: PetscViewerASCIIPrintf(viewer," boomerAMG strength threshold %g\n",jac->as_amg_beta_theta);
1116: if (jac->ams_beta_is_zero_part) {
1117: PetscViewerASCIIPrintf(viewer," compatible subspace projection frequency %d (-1 HYPRE uses default)\n",jac->ams_proj_freq);
1118: }
1119: } else {
1120: PetscViewerASCIIPrintf(viewer," scalar Poisson solver not used (zero-conductivity everywhere) \n");
1121: }
1122: }
1123: return(0);
1124: }
1126: static PetscErrorCode PCSetFromOptions_HYPRE_ADS(PetscOptionItems *PetscOptionsObject,PC pc)
1127: {
1128: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1130: PetscInt n;
1131: PetscBool flag,flag2,flag3,flag4;
1134: PetscOptionsHead(PetscOptionsObject,"HYPRE ADS Options");
1135: PetscOptionsInt("-pc_hypre_ads_print_level","Debugging output level for ADS","None",jac->as_print,&jac->as_print,&flag);
1136: if (flag) PetscStackCallStandard(HYPRE_ADSSetPrintLevel,(jac->hsolver,jac->as_print));
1137: PetscOptionsInt("-pc_hypre_ads_max_iter","Maximum number of ADS multigrid iterations within PCApply","None",jac->as_max_iter,&jac->as_max_iter,&flag);
1138: if (flag) PetscStackCallStandard(HYPRE_ADSSetMaxIter,(jac->hsolver,jac->as_max_iter));
1139: PetscOptionsInt("-pc_hypre_ads_cycle_type","Cycle type for ADS multigrid","None",jac->ads_cycle_type,&jac->ads_cycle_type,&flag);
1140: if (flag) PetscStackCallStandard(HYPRE_ADSSetCycleType,(jac->hsolver,jac->ads_cycle_type));
1141: PetscOptionsReal("-pc_hypre_ads_tol","Error tolerance for ADS multigrid","None",jac->as_tol,&jac->as_tol,&flag);
1142: if (flag) PetscStackCallStandard(HYPRE_ADSSetTol,(jac->hsolver,jac->as_tol));
1143: PetscOptionsInt("-pc_hypre_ads_relax_type","Relaxation type for ADS smoother","None",jac->as_relax_type,&jac->as_relax_type,&flag);
1144: PetscOptionsInt("-pc_hypre_ads_relax_times","Number of relaxation steps for ADS smoother","None",jac->as_relax_times,&jac->as_relax_times,&flag2);
1145: PetscOptionsReal("-pc_hypre_ads_relax_weight","Relaxation weight for ADS smoother","None",jac->as_relax_weight,&jac->as_relax_weight,&flag3);
1146: PetscOptionsReal("-pc_hypre_ads_omega","SSOR coefficient for ADS smoother","None",jac->as_omega,&jac->as_omega,&flag4);
1147: if (flag || flag2 || flag3 || flag4) {
1148: PetscStackCallStandard(HYPRE_ADSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type,
1149: jac->as_relax_times,
1150: jac->as_relax_weight,
1151: jac->as_omega));
1152: }
1153: 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);
1154: n = 5;
1155: PetscOptionsIntArray("-pc_hypre_ads_ams_options","AMG options for AMS solver inside ADS","None",jac->as_amg_alpha_opts,&n,&flag2);
1156: PetscOptionsInt("-pc_hypre_ads_ams_cycle_type","Cycle type for AMS solver inside ADS","None",jac->ams_cycle_type,&jac->ams_cycle_type,&flag3);
1157: if (flag || flag2 || flag3) {
1158: PetscStackCallStandard(HYPRE_ADSSetAMSOptions,(jac->hsolver,jac->ams_cycle_type, /* AMS cycle type */
1159: jac->as_amg_alpha_opts[0], /* AMG coarsen type */
1160: jac->as_amg_alpha_opts[1], /* AMG agg_levels */
1161: jac->as_amg_alpha_opts[2], /* AMG relax_type */
1162: jac->as_amg_alpha_theta,
1163: jac->as_amg_alpha_opts[3], /* AMG interp_type */
1164: jac->as_amg_alpha_opts[4])); /* AMG Pmax */
1165: }
1166: 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);
1167: n = 5;
1168: PetscOptionsIntArray("-pc_hypre_ads_amg_options","AMG options for vector AMG solver inside ADS","None",jac->as_amg_beta_opts,&n,&flag2);
1169: if (flag || flag2) {
1170: PetscStackCallStandard(HYPRE_ADSSetAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */
1171: jac->as_amg_beta_opts[1], /* AMG agg_levels */
1172: jac->as_amg_beta_opts[2], /* AMG relax_type */
1173: jac->as_amg_beta_theta,
1174: jac->as_amg_beta_opts[3], /* AMG interp_type */
1175: jac->as_amg_beta_opts[4])); /* AMG Pmax */
1176: }
1177: PetscOptionsTail();
1178: return(0);
1179: }
1181: static PetscErrorCode PCView_HYPRE_ADS(PC pc,PetscViewer viewer)
1182: {
1183: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1185: PetscBool iascii;
1188: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
1189: if (iascii) {
1190: PetscViewerASCIIPrintf(viewer," HYPRE ADS preconditioning\n");
1191: PetscViewerASCIIPrintf(viewer," subspace iterations per application %d\n",jac->as_max_iter);
1192: PetscViewerASCIIPrintf(viewer," subspace cycle type %d\n",jac->ads_cycle_type);
1193: PetscViewerASCIIPrintf(viewer," subspace iteration tolerance %g\n",jac->as_tol);
1194: PetscViewerASCIIPrintf(viewer," smoother type %d\n",jac->as_relax_type);
1195: PetscViewerASCIIPrintf(viewer," number of smoothing steps %d\n",jac->as_relax_times);
1196: PetscViewerASCIIPrintf(viewer," smoother weight %g\n",jac->as_relax_weight);
1197: PetscViewerASCIIPrintf(viewer," smoother omega %g\n",jac->as_omega);
1198: PetscViewerASCIIPrintf(viewer," AMS solver using boomerAMG\n");
1199: PetscViewerASCIIPrintf(viewer," subspace cycle type %d\n",jac->ams_cycle_type);
1200: PetscViewerASCIIPrintf(viewer," coarsening type %d\n",jac->as_amg_alpha_opts[0]);
1201: PetscViewerASCIIPrintf(viewer," levels of aggressive coarsening %d\n",jac->as_amg_alpha_opts[1]);
1202: PetscViewerASCIIPrintf(viewer," relaxation type %d\n",jac->as_amg_alpha_opts[2]);
1203: PetscViewerASCIIPrintf(viewer," interpolation type %d\n",jac->as_amg_alpha_opts[3]);
1204: PetscViewerASCIIPrintf(viewer," max nonzero elements in interpolation rows %d\n",jac->as_amg_alpha_opts[4]);
1205: PetscViewerASCIIPrintf(viewer," strength threshold %g\n",jac->as_amg_alpha_theta);
1206: PetscViewerASCIIPrintf(viewer," vector Poisson solver using boomerAMG\n");
1207: PetscViewerASCIIPrintf(viewer," coarsening type %d\n",jac->as_amg_beta_opts[0]);
1208: PetscViewerASCIIPrintf(viewer," levels of aggressive coarsening %d\n",jac->as_amg_beta_opts[1]);
1209: PetscViewerASCIIPrintf(viewer," relaxation type %d\n",jac->as_amg_beta_opts[2]);
1210: PetscViewerASCIIPrintf(viewer," interpolation type %d\n",jac->as_amg_beta_opts[3]);
1211: PetscViewerASCIIPrintf(viewer," max nonzero elements in interpolation rows %d\n",jac->as_amg_beta_opts[4]);
1212: PetscViewerASCIIPrintf(viewer," strength threshold %g\n",jac->as_amg_beta_theta);
1213: }
1214: return(0);
1215: }
1217: static PetscErrorCode PCHYPRESetDiscreteGradient_HYPRE(PC pc, Mat G)
1218: {
1219: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1220: PetscBool ishypre;
1224: PetscObjectTypeCompare((PetscObject)G,MATHYPRE,&ishypre);
1225: if (ishypre) {
1226: PetscObjectReference((PetscObject)G);
1227: MatDestroy(&jac->G);
1228: jac->G = G;
1229: } else {
1230: MatDestroy(&jac->G);
1231: MatConvert(G,MATHYPRE,MAT_INITIAL_MATRIX,&jac->G);
1232: }
1233: return(0);
1234: }
1236: /*@
1237: PCHYPRESetDiscreteGradient - Set discrete gradient matrix
1239: Collective on PC
1241: Input Parameters:
1242: + pc - the preconditioning context
1243: - G - the discrete gradient
1245: Level: intermediate
1247: Notes:
1248: G should have as many rows as the number of edges and as many columns as the number of vertices in the mesh
1249: 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
1251: .seealso:
1252: @*/
1253: PetscErrorCode PCHYPRESetDiscreteGradient(PC pc, Mat G)
1254: {
1261: PetscTryMethod(pc,"PCHYPRESetDiscreteGradient_C",(PC,Mat),(pc,G));
1262: return(0);
1263: }
1265: static PetscErrorCode PCHYPRESetDiscreteCurl_HYPRE(PC pc, Mat C)
1266: {
1267: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1268: PetscBool ishypre;
1272: PetscObjectTypeCompare((PetscObject)C,MATHYPRE,&ishypre);
1273: if (ishypre) {
1274: PetscObjectReference((PetscObject)C);
1275: MatDestroy(&jac->C);
1276: jac->C = C;
1277: } else {
1278: MatDestroy(&jac->C);
1279: MatConvert(C,MATHYPRE,MAT_INITIAL_MATRIX,&jac->C);
1280: }
1281: return(0);
1282: }
1284: /*@
1285: PCHYPRESetDiscreteCurl - Set discrete curl matrix
1287: Collective on PC
1289: Input Parameters:
1290: + pc - the preconditioning context
1291: - C - the discrete curl
1293: Level: intermediate
1295: Notes:
1296: C should have as many rows as the number of faces and as many columns as the number of edges in the mesh
1297: 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
1299: .seealso:
1300: @*/
1301: PetscErrorCode PCHYPRESetDiscreteCurl(PC pc, Mat C)
1302: {
1309: PetscTryMethod(pc,"PCHYPRESetDiscreteCurl_C",(PC,Mat),(pc,C));
1310: return(0);
1311: }
1313: static PetscErrorCode PCHYPRESetInterpolations_HYPRE(PC pc, PetscInt dim, Mat RT_PiFull, Mat RT_Pi[], Mat ND_PiFull, Mat ND_Pi[])
1314: {
1315: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1316: PetscBool ishypre;
1318: PetscInt i;
1321: MatDestroy(&jac->RT_PiFull);
1322: MatDestroy(&jac->ND_PiFull);
1323: for (i=0;i<3;++i) {
1324: MatDestroy(&jac->RT_Pi[i]);
1325: MatDestroy(&jac->ND_Pi[i]);
1326: }
1328: jac->dim = dim;
1329: if (RT_PiFull) {
1330: PetscObjectTypeCompare((PetscObject)RT_PiFull,MATHYPRE,&ishypre);
1331: if (ishypre) {
1332: PetscObjectReference((PetscObject)RT_PiFull);
1333: jac->RT_PiFull = RT_PiFull;
1334: } else {
1335: MatConvert(RT_PiFull,MATHYPRE,MAT_INITIAL_MATRIX,&jac->RT_PiFull);
1336: }
1337: }
1338: if (RT_Pi) {
1339: for (i=0;i<dim;++i) {
1340: if (RT_Pi[i]) {
1341: PetscObjectTypeCompare((PetscObject)RT_Pi[i],MATHYPRE,&ishypre);
1342: if (ishypre) {
1343: PetscObjectReference((PetscObject)RT_Pi[i]);
1344: jac->RT_Pi[i] = RT_Pi[i];
1345: } else {
1346: MatConvert(RT_Pi[i],MATHYPRE,MAT_INITIAL_MATRIX,&jac->RT_Pi[i]);
1347: }
1348: }
1349: }
1350: }
1351: if (ND_PiFull) {
1352: PetscObjectTypeCompare((PetscObject)ND_PiFull,MATHYPRE,&ishypre);
1353: if (ishypre) {
1354: PetscObjectReference((PetscObject)ND_PiFull);
1355: jac->ND_PiFull = ND_PiFull;
1356: } else {
1357: MatConvert(ND_PiFull,MATHYPRE,MAT_INITIAL_MATRIX,&jac->ND_PiFull);
1358: }
1359: }
1360: if (ND_Pi) {
1361: for (i=0;i<dim;++i) {
1362: if (ND_Pi[i]) {
1363: PetscObjectTypeCompare((PetscObject)ND_Pi[i],MATHYPRE,&ishypre);
1364: if (ishypre) {
1365: PetscObjectReference((PetscObject)ND_Pi[i]);
1366: jac->ND_Pi[i] = ND_Pi[i];
1367: } else {
1368: MatConvert(ND_Pi[i],MATHYPRE,MAT_INITIAL_MATRIX,&jac->ND_Pi[i]);
1369: }
1370: }
1371: }
1372: }
1374: return(0);
1375: }
1377: /*@
1378: PCHYPRESetInterpolations - Set interpolation matrices for AMS/ADS preconditioner
1380: Collective on PC
1382: Input Parameters:
1383: + pc - the preconditioning context
1384: - dim - the dimension of the problem, only used in AMS
1385: - RT_PiFull - Raviart-Thomas interpolation matrix
1386: - RT_Pi - x/y/z component of Raviart-Thomas interpolation matrix
1387: - ND_PiFull - Nedelec interpolation matrix
1388: - ND_Pi - x/y/z component of Nedelec interpolation matrix
1390: Notes:
1391: For AMS, only Nedelec interpolation matrices are needed, the Raviart-Thomas interpolation matrices can be set to NULL.
1392: For ADS, both type of interpolation matrices are needed.
1393: Level: intermediate
1395: .seealso:
1396: @*/
1397: PetscErrorCode PCHYPRESetInterpolations(PC pc, PetscInt dim, Mat RT_PiFull, Mat RT_Pi[], Mat ND_PiFull, Mat ND_Pi[])
1398: {
1400: PetscInt i;
1404: if (RT_PiFull) {
1407: }
1408: if (RT_Pi) {
1410: for (i=0;i<dim;++i) {
1411: if (RT_Pi[i]) {
1414: }
1415: }
1416: }
1417: if (ND_PiFull) {
1420: }
1421: if (ND_Pi) {
1423: for (i=0;i<dim;++i) {
1424: if (ND_Pi[i]) {
1427: }
1428: }
1429: }
1430: PetscTryMethod(pc,"PCHYPRESetInterpolations_C",(PC,PetscInt,Mat,Mat[],Mat,Mat[]),(pc,dim,RT_PiFull,RT_Pi,ND_PiFull,ND_Pi));
1431: return(0);
1432: }
1434: static PetscErrorCode PCHYPRESetPoissonMatrix_HYPRE(PC pc, Mat A, PetscBool isalpha)
1435: {
1436: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1437: PetscBool ishypre;
1441: PetscObjectTypeCompare((PetscObject)A,MATHYPRE,&ishypre);
1442: if (ishypre) {
1443: if (isalpha) {
1444: PetscObjectReference((PetscObject)A);
1445: MatDestroy(&jac->alpha_Poisson);
1446: jac->alpha_Poisson = A;
1447: } else {
1448: if (A) {
1449: PetscObjectReference((PetscObject)A);
1450: } else {
1451: jac->ams_beta_is_zero = PETSC_TRUE;
1452: }
1453: MatDestroy(&jac->beta_Poisson);
1454: jac->beta_Poisson = A;
1455: }
1456: } else {
1457: if (isalpha) {
1458: MatDestroy(&jac->alpha_Poisson);
1459: MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&jac->alpha_Poisson);
1460: } else {
1461: if (A) {
1462: MatDestroy(&jac->beta_Poisson);
1463: MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&jac->beta_Poisson);
1464: } else {
1465: MatDestroy(&jac->beta_Poisson);
1466: jac->ams_beta_is_zero = PETSC_TRUE;
1467: }
1468: }
1469: }
1470: return(0);
1471: }
1473: /*@
1474: PCHYPRESetAlphaPoissonMatrix - Set vector Poisson matrix
1476: Collective on PC
1478: Input Parameters:
1479: + pc - the preconditioning context
1480: - A - the matrix
1482: Level: intermediate
1484: Notes:
1485: A should be obtained by discretizing the vector valued Poisson problem with linear finite elements
1487: .seealso:
1488: @*/
1489: PetscErrorCode PCHYPRESetAlphaPoissonMatrix(PC pc, Mat A)
1490: {
1497: PetscTryMethod(pc,"PCHYPRESetPoissonMatrix_C",(PC,Mat,PetscBool),(pc,A,PETSC_TRUE));
1498: return(0);
1499: }
1501: /*@
1502: PCHYPRESetBetaPoissonMatrix - Set Poisson matrix
1504: Collective on PC
1506: Input Parameters:
1507: + pc - the preconditioning context
1508: - A - the matrix
1510: Level: intermediate
1512: Notes:
1513: A should be obtained by discretizing the Poisson problem with linear finite elements.
1514: Following HYPRE convention, the scalar Poisson solver of AMS can be turned off by passing NULL.
1516: .seealso:
1517: @*/
1518: PetscErrorCode PCHYPRESetBetaPoissonMatrix(PC pc, Mat A)
1519: {
1524: if (A) {
1527: }
1528: PetscTryMethod(pc,"PCHYPRESetPoissonMatrix_C",(PC,Mat,PetscBool),(pc,A,PETSC_FALSE));
1529: return(0);
1530: }
1532: static PetscErrorCode PCHYPRESetEdgeConstantVectors_HYPRE(PC pc,Vec ozz, Vec zoz, Vec zzo)
1533: {
1534: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1535: PetscErrorCode ierr;
1538: /* throw away any vector if already set */
1539: if (jac->constants[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[0]));
1540: if (jac->constants[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[1]));
1541: if (jac->constants[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[2]));
1542: jac->constants[0] = NULL;
1543: jac->constants[1] = NULL;
1544: jac->constants[2] = NULL;
1545: VecHYPRE_IJVectorCreate(ozz,&jac->constants[0]);
1546: VecHYPRE_IJVectorCopy(ozz,jac->constants[0]);
1547: VecHYPRE_IJVectorCreate(zoz,&jac->constants[1]);
1548: VecHYPRE_IJVectorCopy(zoz,jac->constants[1]);
1549: jac->dim = 2;
1550: if (zzo) {
1551: VecHYPRE_IJVectorCreate(zzo,&jac->constants[2]);
1552: VecHYPRE_IJVectorCopy(zzo,jac->constants[2]);
1553: jac->dim++;
1554: }
1555: return(0);
1556: }
1558: /*@
1559: PCHYPRESetEdgeConstantVectors - Set the representation of the constant vector fields in edge element basis
1561: Collective on PC
1563: Input Parameters:
1564: + pc - the preconditioning context
1565: - ozz - vector representing (1,0,0) (or (1,0) in 2D)
1566: - zoz - vector representing (0,1,0) (or (0,1) in 2D)
1567: - zzo - vector representing (0,0,1) (use NULL in 2D)
1569: Level: intermediate
1571: Notes:
1573: .seealso:
1574: @*/
1575: PetscErrorCode PCHYPRESetEdgeConstantVectors(PC pc, Vec ozz, Vec zoz, Vec zzo)
1576: {
1587: PetscTryMethod(pc,"PCHYPRESetEdgeConstantVectors_C",(PC,Vec,Vec,Vec),(pc,ozz,zoz,zzo));
1588: return(0);
1589: }
1591: static PetscErrorCode PCSetCoordinates_HYPRE(PC pc, PetscInt dim, PetscInt nloc, PetscReal *coords)
1592: {
1593: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1594: Vec tv;
1595: PetscInt i;
1596: PetscErrorCode ierr;
1599: /* throw away any coordinate vector if already set */
1600: if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[0]));
1601: if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[1]));
1602: if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[2]));
1603: jac->dim = dim;
1605: /* compute IJ vector for coordinates */
1606: VecCreate(PetscObjectComm((PetscObject)pc),&tv);
1607: VecSetType(tv,VECSTANDARD);
1608: VecSetSizes(tv,nloc,PETSC_DECIDE);
1609: for (i=0;i<dim;i++) {
1610: PetscScalar *array;
1611: PetscInt j;
1613: VecHYPRE_IJVectorCreate(tv,&jac->coords[i]);
1614: VecGetArray(tv,&array);
1615: for (j=0;j<nloc;j++) {
1616: array[j] = coords[j*dim+i];
1617: }
1618: PetscStackCallStandard(HYPRE_IJVectorSetValues,(jac->coords[i],nloc,NULL,array));
1619: PetscStackCallStandard(HYPRE_IJVectorAssemble,(jac->coords[i]));
1620: VecRestoreArray(tv,&array);
1621: }
1622: VecDestroy(&tv);
1623: return(0);
1624: }
1626: /* ---------------------------------------------------------------------------------*/
1628: static PetscErrorCode PCHYPREGetType_HYPRE(PC pc,const char *name[])
1629: {
1630: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1633: *name = jac->hypre_type;
1634: return(0);
1635: }
1637: static PetscErrorCode PCHYPRESetType_HYPRE(PC pc,const char name[])
1638: {
1639: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
1641: PetscBool flag;
1644: if (jac->hypre_type) {
1645: PetscStrcmp(jac->hypre_type,name,&flag);
1646: if (!flag) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ORDER,"Cannot reset the HYPRE preconditioner type once it has been set");
1647: return(0);
1648: } else {
1649: PetscStrallocpy(name, &jac->hypre_type);
1650: }
1652: jac->maxiter = PETSC_DEFAULT;
1653: jac->tol = PETSC_DEFAULT;
1654: jac->printstatistics = PetscLogPrintInfo;
1656: PetscStrcmp("pilut",jac->hypre_type,&flag);
1657: if (flag) {
1658: MPI_Comm_dup(PetscObjectComm((PetscObject)pc),&(jac->comm_hypre));
1659: PetscStackCallStandard(HYPRE_ParCSRPilutCreate,(jac->comm_hypre,&jac->hsolver));
1660: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Pilut;
1661: pc->ops->view = PCView_HYPRE_Pilut;
1662: jac->destroy = HYPRE_ParCSRPilutDestroy;
1663: jac->setup = HYPRE_ParCSRPilutSetup;
1664: jac->solve = HYPRE_ParCSRPilutSolve;
1665: jac->factorrowsize = PETSC_DEFAULT;
1666: return(0);
1667: }
1668: PetscStrcmp("euclid",jac->hypre_type,&flag);
1669: if (flag) {
1670: MPI_Comm_dup(PetscObjectComm((PetscObject)pc),&(jac->comm_hypre));
1671: PetscStackCallStandard(HYPRE_EuclidCreate,(jac->comm_hypre,&jac->hsolver));
1672: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Euclid;
1673: pc->ops->view = PCView_HYPRE_Euclid;
1674: jac->destroy = HYPRE_EuclidDestroy;
1675: jac->setup = HYPRE_EuclidSetup;
1676: jac->solve = HYPRE_EuclidSolve;
1677: jac->factorrowsize = PETSC_DEFAULT;
1678: jac->eu_level = PETSC_DEFAULT; /* default */
1679: return(0);
1680: }
1681: PetscStrcmp("parasails",jac->hypre_type,&flag);
1682: if (flag) {
1683: MPI_Comm_dup(PetscObjectComm((PetscObject)pc),&(jac->comm_hypre));
1684: PetscStackCallStandard(HYPRE_ParaSailsCreate,(jac->comm_hypre,&jac->hsolver));
1685: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ParaSails;
1686: pc->ops->view = PCView_HYPRE_ParaSails;
1687: jac->destroy = HYPRE_ParaSailsDestroy;
1688: jac->setup = HYPRE_ParaSailsSetup;
1689: jac->solve = HYPRE_ParaSailsSolve;
1690: /* initialize */
1691: jac->nlevels = 1;
1692: jac->threshhold = .1;
1693: jac->filter = .1;
1694: jac->loadbal = 0;
1695: if (PetscLogPrintInfo) jac->logging = (int) PETSC_TRUE;
1696: else jac->logging = (int) PETSC_FALSE;
1698: jac->ruse = (int) PETSC_FALSE;
1699: jac->symt = 0;
1700: PetscStackCallStandard(HYPRE_ParaSailsSetParams,(jac->hsolver,jac->threshhold,jac->nlevels));
1701: PetscStackCallStandard(HYPRE_ParaSailsSetFilter,(jac->hsolver,jac->filter));
1702: PetscStackCallStandard(HYPRE_ParaSailsSetLoadbal,(jac->hsolver,jac->loadbal));
1703: PetscStackCallStandard(HYPRE_ParaSailsSetLogging,(jac->hsolver,jac->logging));
1704: PetscStackCallStandard(HYPRE_ParaSailsSetReuse,(jac->hsolver,jac->ruse));
1705: PetscStackCallStandard(HYPRE_ParaSailsSetSym,(jac->hsolver,jac->symt));
1706: return(0);
1707: }
1708: PetscStrcmp("boomeramg",jac->hypre_type,&flag);
1709: if (flag) {
1710: HYPRE_BoomerAMGCreate(&jac->hsolver);
1711: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_BoomerAMG;
1712: pc->ops->view = PCView_HYPRE_BoomerAMG;
1713: pc->ops->applytranspose = PCApplyTranspose_HYPRE_BoomerAMG;
1714: pc->ops->applyrichardson = PCApplyRichardson_HYPRE_BoomerAMG;
1715: jac->destroy = HYPRE_BoomerAMGDestroy;
1716: jac->setup = HYPRE_BoomerAMGSetup;
1717: jac->solve = HYPRE_BoomerAMGSolve;
1718: jac->applyrichardson = PETSC_FALSE;
1719: /* these defaults match the hypre defaults */
1720: jac->cycletype = 1;
1721: jac->maxlevels = 25;
1722: jac->maxiter = 1;
1723: jac->tol = 0.0; /* tolerance of zero indicates use as preconditioner (suppresses convergence errors) */
1724: jac->truncfactor = 0.0;
1725: jac->strongthreshold = .25;
1726: jac->maxrowsum = .9;
1727: jac->coarsentype = 6;
1728: jac->measuretype = 0;
1729: jac->gridsweeps[0] = jac->gridsweeps[1] = jac->gridsweeps[2] = 1;
1730: jac->smoothtype = -1; /* Not set by default */
1731: jac->smoothnumlevels = 25;
1732: jac->eu_level = 0;
1733: jac->eu_droptolerance = 0;
1734: jac->eu_bj = 0;
1735: jac->relaxtype[0] = jac->relaxtype[1] = 6; /* Defaults to SYMMETRIC since in PETSc we are using a a PC - most likely with CG */
1736: jac->relaxtype[2] = 9; /*G.E. */
1737: jac->relaxweight = 1.0;
1738: jac->outerrelaxweight = 1.0;
1739: jac->relaxorder = 1;
1740: jac->interptype = 0;
1741: jac->agg_nl = 0;
1742: jac->pmax = 0;
1743: jac->truncfactor = 0.0;
1744: jac->agg_num_paths = 1;
1746: jac->nodal_coarsening = 0;
1747: jac->nodal_coarsening_diag = 0;
1748: jac->vec_interp_variant = 0;
1749: jac->vec_interp_qmax = 0;
1750: jac->vec_interp_smooth = PETSC_FALSE;
1751: jac->interp_refine = 0;
1752: jac->nodal_relax = PETSC_FALSE;
1753: jac->nodal_relax_levels = 1;
1754: PetscStackCallStandard(HYPRE_BoomerAMGSetCycleType,(jac->hsolver,jac->cycletype));
1755: PetscStackCallStandard(HYPRE_BoomerAMGSetMaxLevels,(jac->hsolver,jac->maxlevels));
1756: PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,jac->maxiter));
1757: PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,jac->tol));
1758: PetscStackCallStandard(HYPRE_BoomerAMGSetTruncFactor,(jac->hsolver,jac->truncfactor));
1759: PetscStackCallStandard(HYPRE_BoomerAMGSetStrongThreshold,(jac->hsolver,jac->strongthreshold));
1760: PetscStackCallStandard(HYPRE_BoomerAMGSetMaxRowSum,(jac->hsolver,jac->maxrowsum));
1761: PetscStackCallStandard(HYPRE_BoomerAMGSetCoarsenType,(jac->hsolver,jac->coarsentype));
1762: PetscStackCallStandard(HYPRE_BoomerAMGSetMeasureType,(jac->hsolver,jac->measuretype));
1763: PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxOrder,(jac->hsolver, jac->relaxorder));
1764: PetscStackCallStandard(HYPRE_BoomerAMGSetInterpType,(jac->hsolver,jac->interptype));
1765: PetscStackCallStandard(HYPRE_BoomerAMGSetAggNumLevels,(jac->hsolver,jac->agg_nl));
1766: PetscStackCallStandard(HYPRE_BoomerAMGSetPMaxElmts,(jac->hsolver,jac->pmax));
1767: PetscStackCallStandard(HYPRE_BoomerAMGSetNumPaths,(jac->hsolver,jac->agg_num_paths));
1768: PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxType,(jac->hsolver, jac->relaxtype[0])); /*defaults coarse to 9*/
1769: PetscStackCallStandard(HYPRE_BoomerAMGSetNumSweeps,(jac->hsolver, jac->gridsweeps[0])); /*defaults coarse to 1 */
1770: return(0);
1771: }
1772: PetscStrcmp("ams",jac->hypre_type,&flag);
1773: if (flag) {
1774: HYPRE_AMSCreate(&jac->hsolver);
1775: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_AMS;
1776: pc->ops->view = PCView_HYPRE_AMS;
1777: jac->destroy = HYPRE_AMSDestroy;
1778: jac->setup = HYPRE_AMSSetup;
1779: jac->solve = HYPRE_AMSSolve;
1780: jac->coords[0] = NULL;
1781: jac->coords[1] = NULL;
1782: jac->coords[2] = NULL;
1783: /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE AMS */
1784: jac->as_print = 0;
1785: jac->as_max_iter = 1; /* used as a preconditioner */
1786: jac->as_tol = 0.; /* used as a preconditioner */
1787: jac->ams_cycle_type = 13;
1788: /* Smoothing options */
1789: jac->as_relax_type = 2;
1790: jac->as_relax_times = 1;
1791: jac->as_relax_weight = 1.0;
1792: jac->as_omega = 1.0;
1793: /* Vector valued Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
1794: jac->as_amg_alpha_opts[0] = 10;
1795: jac->as_amg_alpha_opts[1] = 1;
1796: jac->as_amg_alpha_opts[2] = 6;
1797: jac->as_amg_alpha_opts[3] = 6;
1798: jac->as_amg_alpha_opts[4] = 4;
1799: jac->as_amg_alpha_theta = 0.25;
1800: /* Scalar Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
1801: jac->as_amg_beta_opts[0] = 10;
1802: jac->as_amg_beta_opts[1] = 1;
1803: jac->as_amg_beta_opts[2] = 6;
1804: jac->as_amg_beta_opts[3] = 6;
1805: jac->as_amg_beta_opts[4] = 4;
1806: jac->as_amg_beta_theta = 0.25;
1807: PetscStackCallStandard(HYPRE_AMSSetPrintLevel,(jac->hsolver,jac->as_print));
1808: PetscStackCallStandard(HYPRE_AMSSetMaxIter,(jac->hsolver,jac->as_max_iter));
1809: PetscStackCallStandard(HYPRE_AMSSetCycleType,(jac->hsolver,jac->ams_cycle_type));
1810: PetscStackCallStandard(HYPRE_AMSSetTol,(jac->hsolver,jac->as_tol));
1811: PetscStackCallStandard(HYPRE_AMSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type,
1812: jac->as_relax_times,
1813: jac->as_relax_weight,
1814: jac->as_omega));
1815: PetscStackCallStandard(HYPRE_AMSSetAlphaAMGOptions,(jac->hsolver,jac->as_amg_alpha_opts[0], /* AMG coarsen type */
1816: jac->as_amg_alpha_opts[1], /* AMG agg_levels */
1817: jac->as_amg_alpha_opts[2], /* AMG relax_type */
1818: jac->as_amg_alpha_theta,
1819: jac->as_amg_alpha_opts[3], /* AMG interp_type */
1820: jac->as_amg_alpha_opts[4])); /* AMG Pmax */
1821: PetscStackCallStandard(HYPRE_AMSSetBetaAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */
1822: jac->as_amg_beta_opts[1], /* AMG agg_levels */
1823: jac->as_amg_beta_opts[2], /* AMG relax_type */
1824: jac->as_amg_beta_theta,
1825: jac->as_amg_beta_opts[3], /* AMG interp_type */
1826: jac->as_amg_beta_opts[4])); /* AMG Pmax */
1827: /* Zero conductivity */
1828: jac->ams_beta_is_zero = PETSC_FALSE;
1829: jac->ams_beta_is_zero_part = PETSC_FALSE;
1830: return(0);
1831: }
1832: PetscStrcmp("ads",jac->hypre_type,&flag);
1833: if (flag) {
1834: HYPRE_ADSCreate(&jac->hsolver);
1835: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ADS;
1836: pc->ops->view = PCView_HYPRE_ADS;
1837: jac->destroy = HYPRE_ADSDestroy;
1838: jac->setup = HYPRE_ADSSetup;
1839: jac->solve = HYPRE_ADSSolve;
1840: jac->coords[0] = NULL;
1841: jac->coords[1] = NULL;
1842: jac->coords[2] = NULL;
1843: /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE ADS */
1844: jac->as_print = 0;
1845: jac->as_max_iter = 1; /* used as a preconditioner */
1846: jac->as_tol = 0.; /* used as a preconditioner */
1847: jac->ads_cycle_type = 13;
1848: /* Smoothing options */
1849: jac->as_relax_type = 2;
1850: jac->as_relax_times = 1;
1851: jac->as_relax_weight = 1.0;
1852: jac->as_omega = 1.0;
1853: /* AMS solver parameters: cycle_type, coarsen type, agg_levels, relax_type, interp_type, Pmax */
1854: jac->ams_cycle_type = 14;
1855: jac->as_amg_alpha_opts[0] = 10;
1856: jac->as_amg_alpha_opts[1] = 1;
1857: jac->as_amg_alpha_opts[2] = 6;
1858: jac->as_amg_alpha_opts[3] = 6;
1859: jac->as_amg_alpha_opts[4] = 4;
1860: jac->as_amg_alpha_theta = 0.25;
1861: /* Vector Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
1862: jac->as_amg_beta_opts[0] = 10;
1863: jac->as_amg_beta_opts[1] = 1;
1864: jac->as_amg_beta_opts[2] = 6;
1865: jac->as_amg_beta_opts[3] = 6;
1866: jac->as_amg_beta_opts[4] = 4;
1867: jac->as_amg_beta_theta = 0.25;
1868: PetscStackCallStandard(HYPRE_ADSSetPrintLevel,(jac->hsolver,jac->as_print));
1869: PetscStackCallStandard(HYPRE_ADSSetMaxIter,(jac->hsolver,jac->as_max_iter));
1870: PetscStackCallStandard(HYPRE_ADSSetCycleType,(jac->hsolver,jac->ams_cycle_type));
1871: PetscStackCallStandard(HYPRE_ADSSetTol,(jac->hsolver,jac->as_tol));
1872: PetscStackCallStandard(HYPRE_ADSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type,
1873: jac->as_relax_times,
1874: jac->as_relax_weight,
1875: jac->as_omega));
1876: PetscStackCallStandard(HYPRE_ADSSetAMSOptions,(jac->hsolver,jac->ams_cycle_type, /* AMG coarsen type */
1877: jac->as_amg_alpha_opts[0], /* AMG coarsen type */
1878: jac->as_amg_alpha_opts[1], /* AMG agg_levels */
1879: jac->as_amg_alpha_opts[2], /* AMG relax_type */
1880: jac->as_amg_alpha_theta,
1881: jac->as_amg_alpha_opts[3], /* AMG interp_type */
1882: jac->as_amg_alpha_opts[4])); /* AMG Pmax */
1883: PetscStackCallStandard(HYPRE_ADSSetAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */
1884: jac->as_amg_beta_opts[1], /* AMG agg_levels */
1885: jac->as_amg_beta_opts[2], /* AMG relax_type */
1886: jac->as_amg_beta_theta,
1887: jac->as_amg_beta_opts[3], /* AMG interp_type */
1888: jac->as_amg_beta_opts[4])); /* AMG Pmax */
1889: return(0);
1890: }
1891: PetscFree(jac->hypre_type);
1893: jac->hypre_type = NULL;
1894: SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown HYPRE preconditioner %s; Choices are euclid, pilut, parasails, boomeramg, ams",name);
1895: return(0);
1896: }
1898: /*
1899: It only gets here if the HYPRE type has not been set before the call to
1900: ...SetFromOptions() which actually is most of the time
1901: */
1902: static PetscErrorCode PCSetFromOptions_HYPRE(PetscOptionItems *PetscOptionsObject,PC pc)
1903: {
1905: PetscInt indx;
1906: const char *type[] = {"euclid","pilut","parasails","boomeramg","ams","ads"};
1907: PetscBool flg;
1910: PetscOptionsHead(PetscOptionsObject,"HYPRE preconditioner options");
1911: PetscOptionsEList("-pc_hypre_type","HYPRE preconditioner type","PCHYPRESetType",type,4,"boomeramg",&indx,&flg);
1912: if (flg) {
1913: PCHYPRESetType_HYPRE(pc,type[indx]);
1914: } else {
1915: PCHYPRESetType_HYPRE(pc,"boomeramg");
1916: }
1917: if (pc->ops->setfromoptions) {
1918: pc->ops->setfromoptions(PetscOptionsObject,pc);
1919: }
1920: PetscOptionsTail();
1921: return(0);
1922: }
1924: /*@C
1925: PCHYPRESetType - Sets which hypre preconditioner you wish to use
1927: Input Parameters:
1928: + pc - the preconditioner context
1929: - name - either euclid, pilut, parasails, boomeramg, ams, ads
1931: Options Database Keys:
1932: -pc_hypre_type - One of euclid, pilut, parasails, boomeramg, ams, ads
1934: Level: intermediate
1936: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
1937: PCHYPRE
1939: @*/
1940: PetscErrorCode PCHYPRESetType(PC pc,const char name[])
1941: {
1947: PetscTryMethod(pc,"PCHYPRESetType_C",(PC,const char[]),(pc,name));
1948: return(0);
1949: }
1951: /*@C
1952: PCHYPREGetType - Gets which hypre preconditioner you are using
1954: Input Parameter:
1955: . pc - the preconditioner context
1957: Output Parameter:
1958: . name - either euclid, pilut, parasails, boomeramg, ams, ads
1960: Level: intermediate
1962: .seealso: PCCreate(), PCHYPRESetType(), PCType (for list of available types), PC,
1963: PCHYPRE
1965: @*/
1966: PetscErrorCode PCHYPREGetType(PC pc,const char *name[])
1967: {
1973: PetscTryMethod(pc,"PCHYPREGetType_C",(PC,const char*[]),(pc,name));
1974: return(0);
1975: }
1977: /*MC
1978: PCHYPRE - Allows you to use the matrix element based preconditioners in the LLNL package hypre
1980: Options Database Keys:
1981: + -pc_hypre_type - One of euclid, pilut, parasails, boomeramg, ams, ads
1982: - Too many others to list, run with -pc_type hypre -pc_hypre_type XXX -help to see options for the XXX
1983: preconditioner
1985: Level: intermediate
1987: Notes:
1988: Apart from pc_hypre_type (for which there is PCHYPRESetType()),
1989: the many hypre options can ONLY be set via the options database (e.g. the command line
1990: or with PetscOptionsSetValue(), there are no functions to set them)
1992: The options -pc_hypre_boomeramg_max_iter and -pc_hypre_boomeramg_tol refer to the number of iterations
1993: (V-cycles) and tolerance that boomeramg does EACH time it is called. So for example, if
1994: -pc_hypre_boomeramg_max_iter is set to 2 then 2-V-cycles are being used to define the preconditioner
1995: (-pc_hypre_boomeramg_tol should be set to 0.0 - the default - to strictly use a fixed number of
1996: iterations per hypre call). -ksp_max_it and -ksp_rtol STILL determine the total number of iterations
1997: and tolerance for the Krylov solver. For example, if -pc_hypre_boomeramg_max_iter is 2 and -ksp_max_it is 10
1998: then AT MOST twenty V-cycles of boomeramg will be called.
2000: Note that the option -pc_hypre_boomeramg_relax_type_all defaults to symmetric relaxation
2001: (symmetric-SOR/Jacobi), which is required for Krylov solvers like CG that expect symmetry.
2002: Otherwise, you may want to use -pc_hypre_boomeramg_relax_type_all SOR/Jacobi.
2003: If you wish to use BoomerAMG WITHOUT a Krylov method use -ksp_type richardson NOT -ksp_type preonly
2004: and use -ksp_max_it to control the number of V-cycles.
2005: (see the PETSc FAQ.html at the PETSc website under the Documentation tab).
2007: 2007-02-03 Using HYPRE-1.11.1b, the routine HYPRE_BoomerAMGSolveT and the option
2008: -pc_hypre_parasails_reuse were failing with SIGSEGV. Dalcin L.
2010: MatSetNearNullSpace() - if you provide a near null space to your matrix it is ignored by hypre UNLESS you also use
2011: the two options:
2012: + -pc_hypre_boomeramg_nodal_coarsen <n> - where n is from 1 to 6 (see HYPRE_BOOMERAMGSetNodal())
2013: - -pc_hypre_boomeramg_vec_interp_variant <v> where v is from 1 to 3 (see HYPRE_BoomerAMGSetInterpVecVariant())
2015: Depending on the linear system you may see the same or different convergence depending on the values you use.
2017: See PCPFMG for access to the hypre Struct PFMG solver
2019: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
2020: PCHYPRESetType(), PCPFMG
2022: M*/
2024: PETSC_EXTERN PetscErrorCode PCCreate_HYPRE(PC pc)
2025: {
2026: PC_HYPRE *jac;
2030: PetscNewLog(pc,&jac);
2032: pc->data = jac;
2033: pc->ops->reset = PCReset_HYPRE;
2034: pc->ops->destroy = PCDestroy_HYPRE;
2035: pc->ops->setfromoptions = PCSetFromOptions_HYPRE;
2036: pc->ops->setup = PCSetUp_HYPRE;
2037: pc->ops->apply = PCApply_HYPRE;
2038: jac->comm_hypre = MPI_COMM_NULL;
2039: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetType_C",PCHYPRESetType_HYPRE);
2040: PetscObjectComposeFunction((PetscObject)pc,"PCHYPREGetType_C",PCHYPREGetType_HYPRE);
2041: PetscObjectComposeFunction((PetscObject)pc,"PCSetCoordinates_C",PCSetCoordinates_HYPRE);
2042: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteGradient_C",PCHYPRESetDiscreteGradient_HYPRE);
2043: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteCurl_C",PCHYPRESetDiscreteCurl_HYPRE);
2044: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetInterpolations_C",PCHYPRESetInterpolations_HYPRE);
2045: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetEdgeConstantVectors_C",PCHYPRESetEdgeConstantVectors_HYPRE);
2046: PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetPoissonMatrix_C",PCHYPRESetPoissonMatrix_HYPRE);
2047: return(0);
2048: }
2050: /* ---------------------------------------------------------------------------------------------------------------------------------*/
2052: typedef struct {
2053: MPI_Comm hcomm; /* does not share comm with HYPRE_StructMatrix because need to create solver before getting matrix */
2054: HYPRE_StructSolver hsolver;
2056: /* keep copy of PFMG options used so may view them */
2057: PetscInt its;
2058: double tol;
2059: PetscInt relax_type;
2060: PetscInt rap_type;
2061: PetscInt num_pre_relax,num_post_relax;
2062: PetscInt max_levels;
2063: } PC_PFMG;
2065: PetscErrorCode PCDestroy_PFMG(PC pc)
2066: {
2068: PC_PFMG *ex = (PC_PFMG*) pc->data;
2071: if (ex->hsolver) PetscStackCallStandard(HYPRE_StructPFMGDestroy,(ex->hsolver));
2072: MPI_Comm_free(&ex->hcomm);
2073: PetscFree(pc->data);
2074: return(0);
2075: }
2077: static const char *PFMGRelaxType[] = {"Jacobi","Weighted-Jacobi","symmetric-Red/Black-Gauss-Seidel","Red/Black-Gauss-Seidel"};
2078: static const char *PFMGRAPType[] = {"Galerkin","non-Galerkin"};
2080: PetscErrorCode PCView_PFMG(PC pc,PetscViewer viewer)
2081: {
2083: PetscBool iascii;
2084: PC_PFMG *ex = (PC_PFMG*) pc->data;
2087: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
2088: if (iascii) {
2089: PetscViewerASCIIPrintf(viewer," HYPRE PFMG preconditioning\n");
2090: PetscViewerASCIIPrintf(viewer," max iterations %d\n",ex->its);
2091: PetscViewerASCIIPrintf(viewer," tolerance %g\n",ex->tol);
2092: PetscViewerASCIIPrintf(viewer," relax type %s\n",PFMGRelaxType[ex->relax_type]);
2093: PetscViewerASCIIPrintf(viewer," RAP type %s\n",PFMGRAPType[ex->rap_type]);
2094: PetscViewerASCIIPrintf(viewer," number pre-relax %d post-relax %d\n",ex->num_pre_relax,ex->num_post_relax);
2095: PetscViewerASCIIPrintf(viewer," max levels %d\n",ex->max_levels);
2096: }
2097: return(0);
2098: }
2100: PetscErrorCode PCSetFromOptions_PFMG(PetscOptionItems *PetscOptionsObject,PC pc)
2101: {
2103: PC_PFMG *ex = (PC_PFMG*) pc->data;
2104: PetscBool flg = PETSC_FALSE;
2107: PetscOptionsHead(PetscOptionsObject,"PFMG options");
2108: PetscOptionsBool("-pc_pfmg_print_statistics","Print statistics","HYPRE_StructPFMGSetPrintLevel",flg,&flg,NULL);
2109: if (flg) {
2110: int level=3;
2111: PetscStackCallStandard(HYPRE_StructPFMGSetPrintLevel,(ex->hsolver,level));
2112: }
2113: PetscOptionsInt("-pc_pfmg_its","Number of iterations of PFMG to use as preconditioner","HYPRE_StructPFMGSetMaxIter",ex->its,&ex->its,NULL);
2114: PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,(ex->hsolver,ex->its));
2115: PetscOptionsInt("-pc_pfmg_num_pre_relax","Number of smoothing steps before coarse grid","HYPRE_StructPFMGSetNumPreRelax",ex->num_pre_relax,&ex->num_pre_relax,NULL);
2116: PetscStackCallStandard(HYPRE_StructPFMGSetNumPreRelax,(ex->hsolver,ex->num_pre_relax));
2117: PetscOptionsInt("-pc_pfmg_num_post_relax","Number of smoothing steps after coarse grid","HYPRE_StructPFMGSetNumPostRelax",ex->num_post_relax,&ex->num_post_relax,NULL);
2118: PetscStackCallStandard(HYPRE_StructPFMGSetNumPostRelax,(ex->hsolver,ex->num_post_relax));
2120: PetscOptionsInt("-pc_pfmg_max_levels","Max Levels for MG hierarchy","HYPRE_StructPFMGSetMaxLevels",ex->max_levels,&ex->max_levels,NULL);
2121: PetscStackCallStandard(HYPRE_StructPFMGSetMaxLevels,(ex->hsolver,ex->max_levels));
2123: PetscOptionsReal("-pc_pfmg_tol","Tolerance of PFMG","HYPRE_StructPFMGSetTol",ex->tol,&ex->tol,NULL);
2124: PetscStackCallStandard(HYPRE_StructPFMGSetTol,(ex->hsolver,ex->tol));
2125: 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);
2126: PetscStackCallStandard(HYPRE_StructPFMGSetRelaxType,(ex->hsolver, ex->relax_type));
2127: PetscOptionsEList("-pc_pfmg_rap_type","RAP type","HYPRE_StructPFMGSetRAPType",PFMGRAPType,ALEN(PFMGRAPType),PFMGRAPType[ex->rap_type],&ex->rap_type,NULL);
2128: PetscStackCallStandard(HYPRE_StructPFMGSetRAPType,(ex->hsolver, ex->rap_type));
2129: PetscOptionsTail();
2130: return(0);
2131: }
2133: PetscErrorCode PCApply_PFMG(PC pc,Vec x,Vec y)
2134: {
2135: PetscErrorCode ierr;
2136: PC_PFMG *ex = (PC_PFMG*) pc->data;
2137: PetscScalar *yy;
2138: const PetscScalar *xx;
2139: PetscInt ilower[3],iupper[3];
2140: Mat_HYPREStruct *mx = (Mat_HYPREStruct*)(pc->pmat->data);
2143: PetscCitationsRegister(hypreCitation,&cite);
2144: DMDAGetCorners(mx->da,&ilower[0],&ilower[1],&ilower[2],&iupper[0],&iupper[1],&iupper[2]);
2145: iupper[0] += ilower[0] - 1;
2146: iupper[1] += ilower[1] - 1;
2147: iupper[2] += ilower[2] - 1;
2149: /* copy x values over to hypre */
2150: PetscStackCallStandard(HYPRE_StructVectorSetConstantValues,(mx->hb,0.0));
2151: VecGetArrayRead(x,&xx);
2152: PetscStackCallStandard(HYPRE_StructVectorSetBoxValues,(mx->hb,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,(PetscScalar*)xx));
2153: VecRestoreArrayRead(x,&xx);
2154: PetscStackCallStandard(HYPRE_StructVectorAssemble,(mx->hb));
2155: PetscStackCallStandard(HYPRE_StructPFMGSolve,(ex->hsolver,mx->hmat,mx->hb,mx->hx));
2157: /* copy solution values back to PETSc */
2158: VecGetArray(y,&yy);
2159: PetscStackCallStandard(HYPRE_StructVectorGetBoxValues,(mx->hx,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,yy));
2160: VecRestoreArray(y,&yy);
2161: return(0);
2162: }
2164: 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)
2165: {
2166: PC_PFMG *jac = (PC_PFMG*)pc->data;
2168: PetscInt oits;
2171: PetscCitationsRegister(hypreCitation,&cite);
2172: PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,(jac->hsolver,its*jac->its));
2173: PetscStackCallStandard(HYPRE_StructPFMGSetTol,(jac->hsolver,rtol));
2175: PCApply_PFMG(pc,b,y);
2176: PetscStackCallStandard(HYPRE_StructPFMGGetNumIterations,(jac->hsolver,(HYPRE_Int *)&oits));
2177: *outits = oits;
2178: if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
2179: else *reason = PCRICHARDSON_CONVERGED_RTOL;
2180: PetscStackCallStandard(HYPRE_StructPFMGSetTol,(jac->hsolver,jac->tol));
2181: PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,(jac->hsolver,jac->its));
2182: return(0);
2183: }
2186: PetscErrorCode PCSetUp_PFMG(PC pc)
2187: {
2188: PetscErrorCode ierr;
2189: PC_PFMG *ex = (PC_PFMG*) pc->data;
2190: Mat_HYPREStruct *mx = (Mat_HYPREStruct*)(pc->pmat->data);
2191: PetscBool flg;
2194: PetscObjectTypeCompare((PetscObject)pc->pmat,MATHYPRESTRUCT,&flg);
2195: if (!flg) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Must use MATHYPRESTRUCT with this preconditioner");
2197: /* create the hypre solver object and set its information */
2198: if (ex->hsolver) PetscStackCallStandard(HYPRE_StructPFMGDestroy,(ex->hsolver));
2199: PetscStackCallStandard(HYPRE_StructPFMGCreate,(ex->hcomm,&ex->hsolver));
2200: PetscStackCallStandard(HYPRE_StructPFMGSetup,(ex->hsolver,mx->hmat,mx->hb,mx->hx));
2201: PetscStackCallStandard(HYPRE_StructPFMGSetZeroGuess,(ex->hsolver));
2202: return(0);
2203: }
2205: /*MC
2206: PCPFMG - the hypre PFMG multigrid solver
2208: Level: advanced
2210: Options Database:
2211: + -pc_pfmg_its <its> number of iterations of PFMG to use as preconditioner
2212: . -pc_pfmg_num_pre_relax <steps> number of smoothing steps before coarse grid
2213: . -pc_pfmg_num_post_relax <steps> number of smoothing steps after coarse grid
2214: . -pc_pfmg_tol <tol> tolerance of PFMG
2215: . -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
2216: - -pc_pfmg_rap_type - type of coarse matrix generation, one of Galerkin,non-Galerkin
2218: Notes:
2219: This is for CELL-centered descretizations
2221: This must be used with the MATHYPRESTRUCT matrix type.
2222: This is less general than in hypre, it supports only one block per process defined by a PETSc DMDA.
2224: .seealso: PCMG, MATHYPRESTRUCT
2225: M*/
2227: PETSC_EXTERN PetscErrorCode PCCreate_PFMG(PC pc)
2228: {
2230: PC_PFMG *ex;
2233: PetscNew(&ex); \
2234: pc->data = ex;
2236: ex->its = 1;
2237: ex->tol = 1.e-8;
2238: ex->relax_type = 1;
2239: ex->rap_type = 0;
2240: ex->num_pre_relax = 1;
2241: ex->num_post_relax = 1;
2242: ex->max_levels = 0;
2244: pc->ops->setfromoptions = PCSetFromOptions_PFMG;
2245: pc->ops->view = PCView_PFMG;
2246: pc->ops->destroy = PCDestroy_PFMG;
2247: pc->ops->apply = PCApply_PFMG;
2248: pc->ops->applyrichardson = PCApplyRichardson_PFMG;
2249: pc->ops->setup = PCSetUp_PFMG;
2251: MPI_Comm_dup(PetscObjectComm((PetscObject)pc),&(ex->hcomm));
2252: PetscStackCallStandard(HYPRE_StructPFMGCreate,(ex->hcomm,&ex->hsolver));
2253: return(0);
2254: }
2256: /* ---------------------------------------------------------------------------------------------------------------------------------------------------*/
2258: /* we know we are working with a HYPRE_SStructMatrix */
2259: typedef struct {
2260: MPI_Comm hcomm; /* does not share comm with HYPRE_SStructMatrix because need to create solver before getting matrix */
2261: HYPRE_SStructSolver ss_solver;
2263: /* keep copy of SYSPFMG options used so may view them */
2264: PetscInt its;
2265: double tol;
2266: PetscInt relax_type;
2267: PetscInt num_pre_relax,num_post_relax;
2268: } PC_SysPFMG;
2270: PetscErrorCode PCDestroy_SysPFMG(PC pc)
2271: {
2273: PC_SysPFMG *ex = (PC_SysPFMG*) pc->data;
2276: if (ex->ss_solver) PetscStackCallStandard(HYPRE_SStructSysPFMGDestroy,(ex->ss_solver));
2277: MPI_Comm_free(&ex->hcomm);
2278: PetscFree(pc->data);
2279: return(0);
2280: }
2282: static const char *SysPFMGRelaxType[] = {"Weighted-Jacobi","Red/Black-Gauss-Seidel"};
2284: PetscErrorCode PCView_SysPFMG(PC pc,PetscViewer viewer)
2285: {
2287: PetscBool iascii;
2288: PC_SysPFMG *ex = (PC_SysPFMG*) pc->data;
2291: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
2292: if (iascii) {
2293: PetscViewerASCIIPrintf(viewer," HYPRE SysPFMG preconditioning\n");
2294: PetscViewerASCIIPrintf(viewer," max iterations %d\n",ex->its);
2295: PetscViewerASCIIPrintf(viewer," tolerance %g\n",ex->tol);
2296: PetscViewerASCIIPrintf(viewer," relax type %s\n",PFMGRelaxType[ex->relax_type]);
2297: PetscViewerASCIIPrintf(viewer," number pre-relax %d post-relax %d\n",ex->num_pre_relax,ex->num_post_relax);
2298: }
2299: return(0);
2300: }
2302: PetscErrorCode PCSetFromOptions_SysPFMG(PetscOptionItems *PetscOptionsObject,PC pc)
2303: {
2305: PC_SysPFMG *ex = (PC_SysPFMG*) pc->data;
2306: PetscBool flg = PETSC_FALSE;
2309: PetscOptionsHead(PetscOptionsObject,"SysPFMG options");
2310: PetscOptionsBool("-pc_syspfmg_print_statistics","Print statistics","HYPRE_SStructSysPFMGSetPrintLevel",flg,&flg,NULL);
2311: if (flg) {
2312: int level=3;
2313: PetscStackCallStandard(HYPRE_SStructSysPFMGSetPrintLevel,(ex->ss_solver,level));
2314: }
2315: PetscOptionsInt("-pc_syspfmg_its","Number of iterations of SysPFMG to use as preconditioner","HYPRE_SStructSysPFMGSetMaxIter",ex->its,&ex->its,NULL);
2316: PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,(ex->ss_solver,ex->its));
2317: PetscOptionsInt("-pc_syspfmg_num_pre_relax","Number of smoothing steps before coarse grid","HYPRE_SStructSysPFMGSetNumPreRelax",ex->num_pre_relax,&ex->num_pre_relax,NULL);
2318: PetscStackCallStandard(HYPRE_SStructSysPFMGSetNumPreRelax,(ex->ss_solver,ex->num_pre_relax));
2319: PetscOptionsInt("-pc_syspfmg_num_post_relax","Number of smoothing steps after coarse grid","HYPRE_SStructSysPFMGSetNumPostRelax",ex->num_post_relax,&ex->num_post_relax,NULL);
2320: PetscStackCallStandard(HYPRE_SStructSysPFMGSetNumPostRelax,(ex->ss_solver,ex->num_post_relax));
2322: PetscOptionsReal("-pc_syspfmg_tol","Tolerance of SysPFMG","HYPRE_SStructSysPFMGSetTol",ex->tol,&ex->tol,NULL);
2323: PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,(ex->ss_solver,ex->tol));
2324: PetscOptionsEList("-pc_syspfmg_relax_type","Relax type for the up and down cycles","HYPRE_SStructSysPFMGSetRelaxType",SysPFMGRelaxType,ALEN(SysPFMGRelaxType),SysPFMGRelaxType[ex->relax_type],&ex->relax_type,NULL);
2325: PetscStackCallStandard(HYPRE_SStructSysPFMGSetRelaxType,(ex->ss_solver, ex->relax_type));
2326: PetscOptionsTail();
2327: return(0);
2328: }
2330: PetscErrorCode PCApply_SysPFMG(PC pc,Vec x,Vec y)
2331: {
2332: PetscErrorCode ierr;
2333: PC_SysPFMG *ex = (PC_SysPFMG*) pc->data;
2334: PetscScalar *yy;
2335: const PetscScalar *xx;
2336: PetscInt ilower[3],iupper[3];
2337: Mat_HYPRESStruct *mx = (Mat_HYPRESStruct*)(pc->pmat->data);
2338: PetscInt ordering= mx->dofs_order;
2339: PetscInt nvars = mx->nvars;
2340: PetscInt part = 0;
2341: PetscInt size;
2342: PetscInt i;
2345: PetscCitationsRegister(hypreCitation,&cite);
2346: DMDAGetCorners(mx->da,&ilower[0],&ilower[1],&ilower[2],&iupper[0],&iupper[1],&iupper[2]);
2347: iupper[0] += ilower[0] - 1;
2348: iupper[1] += ilower[1] - 1;
2349: iupper[2] += ilower[2] - 1;
2351: size = 1;
2352: for (i= 0; i< 3; i++) size *= (iupper[i]-ilower[i]+1);
2354: /* copy x values over to hypre for variable ordering */
2355: if (ordering) {
2356: PetscStackCallStandard(HYPRE_SStructVectorSetConstantValues,(mx->ss_b,0.0));
2357: VecGetArrayRead(x,&xx);
2358: for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorSetBoxValues,(mx->ss_b,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,(PetscScalar*)xx+(size*i)));
2359: VecRestoreArrayRead(x,&xx);
2360: PetscStackCallStandard(HYPRE_SStructVectorAssemble,(mx->ss_b));
2361: PetscStackCallStandard(HYPRE_SStructMatrixMatvec,(1.0,mx->ss_mat,mx->ss_b,0.0,mx->ss_x));
2362: PetscStackCallStandard(HYPRE_SStructSysPFMGSolve,(ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x));
2364: /* copy solution values back to PETSc */
2365: VecGetArray(y,&yy);
2366: for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorGetBoxValues,(mx->ss_x,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,yy+(size*i)));
2367: VecRestoreArray(y,&yy);
2368: } else { /* nodal ordering must be mapped to variable ordering for sys_pfmg */
2369: PetscScalar *z;
2370: PetscInt j, k;
2372: PetscMalloc1(nvars*size,&z);
2373: PetscStackCallStandard(HYPRE_SStructVectorSetConstantValues,(mx->ss_b,0.0));
2374: VecGetArrayRead(x,&xx);
2376: /* transform nodal to hypre's variable ordering for sys_pfmg */
2377: for (i= 0; i< size; i++) {
2378: k= i*nvars;
2379: for (j= 0; j< nvars; j++) z[j*size+i]= xx[k+j];
2380: }
2381: for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorSetBoxValues,(mx->ss_b,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,z+(size*i)));
2382: VecRestoreArrayRead(x,&xx);
2383: PetscStackCallStandard(HYPRE_SStructVectorAssemble,(mx->ss_b));
2384: PetscStackCallStandard(HYPRE_SStructSysPFMGSolve,(ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x));
2386: /* copy solution values back to PETSc */
2387: VecGetArray(y,&yy);
2388: for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorGetBoxValues,(mx->ss_x,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,z+(size*i)));
2389: /* transform hypre's variable ordering for sys_pfmg to nodal ordering */
2390: for (i= 0; i< size; i++) {
2391: k= i*nvars;
2392: for (j= 0; j< nvars; j++) yy[k+j]= z[j*size+i];
2393: }
2394: VecRestoreArray(y,&yy);
2395: PetscFree(z);
2396: }
2397: return(0);
2398: }
2400: 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)
2401: {
2402: PC_SysPFMG *jac = (PC_SysPFMG*)pc->data;
2404: PetscInt oits;
2407: PetscCitationsRegister(hypreCitation,&cite);
2408: PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,(jac->ss_solver,its*jac->its));
2409: PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,(jac->ss_solver,rtol));
2410: PCApply_SysPFMG(pc,b,y);
2411: PetscStackCallStandard(HYPRE_SStructSysPFMGGetNumIterations,(jac->ss_solver,(HYPRE_Int *)&oits));
2412: *outits = oits;
2413: if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
2414: else *reason = PCRICHARDSON_CONVERGED_RTOL;
2415: PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,(jac->ss_solver,jac->tol));
2416: PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,(jac->ss_solver,jac->its));
2417: return(0);
2418: }
2420: PetscErrorCode PCSetUp_SysPFMG(PC pc)
2421: {
2422: PetscErrorCode ierr;
2423: PC_SysPFMG *ex = (PC_SysPFMG*) pc->data;
2424: Mat_HYPRESStruct *mx = (Mat_HYPRESStruct*)(pc->pmat->data);
2425: PetscBool flg;
2428: PetscObjectTypeCompare((PetscObject)pc->pmat,MATHYPRESSTRUCT,&flg);
2429: if (!flg) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Must use MATHYPRESSTRUCT with this preconditioner");
2431: /* create the hypre sstruct solver object and set its information */
2432: if (ex->ss_solver) PetscStackCallStandard(HYPRE_SStructSysPFMGDestroy,(ex->ss_solver));
2433: PetscStackCallStandard(HYPRE_SStructSysPFMGCreate,(ex->hcomm,&ex->ss_solver));
2434: PetscStackCallStandard(HYPRE_SStructSysPFMGSetZeroGuess,(ex->ss_solver));
2435: PetscStackCallStandard(HYPRE_SStructSysPFMGSetup,(ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x));
2436: return(0);
2437: }
2439: /*MC
2440: PCSysPFMG - the hypre SysPFMG multigrid solver
2442: Level: advanced
2444: Options Database:
2445: + -pc_syspfmg_its <its> number of iterations of SysPFMG to use as preconditioner
2446: . -pc_syspfmg_num_pre_relax <steps> number of smoothing steps before coarse grid
2447: . -pc_syspfmg_num_post_relax <steps> number of smoothing steps after coarse grid
2448: . -pc_syspfmg_tol <tol> tolerance of SysPFMG
2449: . -pc_syspfmg_relax_type -relaxation type for the up and down cycles, one of Weighted-Jacobi,Red/Black-Gauss-Seidel
2451: Notes:
2452: This is for CELL-centered descretizations
2454: This must be used with the MATHYPRESSTRUCT matrix type.
2455: This is less general than in hypre, it supports only one part, and one block per process defined by a PETSc DMDA.
2456: Also, only cell-centered variables.
2458: .seealso: PCMG, MATHYPRESSTRUCT
2459: M*/
2461: PETSC_EXTERN PetscErrorCode PCCreate_SysPFMG(PC pc)
2462: {
2464: PC_SysPFMG *ex;
2467: PetscNew(&ex); \
2468: pc->data = ex;
2470: ex->its = 1;
2471: ex->tol = 1.e-8;
2472: ex->relax_type = 1;
2473: ex->num_pre_relax = 1;
2474: ex->num_post_relax = 1;
2476: pc->ops->setfromoptions = PCSetFromOptions_SysPFMG;
2477: pc->ops->view = PCView_SysPFMG;
2478: pc->ops->destroy = PCDestroy_SysPFMG;
2479: pc->ops->apply = PCApply_SysPFMG;
2480: pc->ops->applyrichardson = PCApplyRichardson_SysPFMG;
2481: pc->ops->setup = PCSetUp_SysPFMG;
2483: MPI_Comm_dup(PetscObjectComm((PetscObject)pc),&(ex->hcomm));
2484: PetscStackCallStandard(HYPRE_SStructSysPFMGCreate,(ex->hcomm,&ex->ss_solver));
2485: return(0);
2486: }