Actual source code: hypre.c
1: /*
2: Provides an interface to the LLNL package hypre
3: */
5: #include <petscpkg_version.h>
6: #include <petsc/private/pcimpl.h>
7: /* this include is needed ONLY to allow access to the private data inside the Mat object specific to hypre */
8: #include <petsc/private/matimpl.h>
9: #include <petsc/private/vecimpl.h>
10: #include <../src/vec/vec/impls/hypre/vhyp.h>
11: #include <../src/mat/impls/hypre/mhypre.h>
12: #include <../src/dm/impls/da/hypre/mhyp.h>
13: #include <_hypre_parcsr_ls.h>
14: #include <petscmathypre.h>
16: #if defined(PETSC_HAVE_HYPRE_DEVICE)
17: #include <petsc/private/deviceimpl.h>
18: #endif
20: static PetscBool cite = PETSC_FALSE;
21: static const char hypreCitation[] = "@manual{hypre-web-page,\n title = {{\\sl hypre}: High Performance Preconditioners},\n organization = {Lawrence Livermore National Laboratory},\n note = "
22: "{\\url{https://www.llnl.gov/casc/hypre}}\n}\n";
24: /*
25: Private context (data structure) for the preconditioner.
26: */
27: typedef struct {
28: HYPRE_Solver hsolver;
29: Mat hpmat; /* MatHYPRE */
31: HYPRE_Int (*destroy)(HYPRE_Solver);
32: HYPRE_Int (*solve)(HYPRE_Solver, HYPRE_ParCSRMatrix, HYPRE_ParVector, HYPRE_ParVector);
33: HYPRE_Int (*setup)(HYPRE_Solver, HYPRE_ParCSRMatrix, HYPRE_ParVector, HYPRE_ParVector);
35: MPI_Comm comm_hypre;
36: char *hypre_type;
38: /* options for Pilut and BoomerAMG*/
39: PetscInt maxiter;
40: PetscReal tol;
42: /* options for Pilut */
43: PetscInt factorrowsize;
45: /* options for ParaSails */
46: PetscInt nlevels;
47: PetscReal threshold;
48: PetscReal filter;
49: PetscReal loadbal;
50: PetscInt logging;
51: PetscInt ruse;
52: PetscInt symt;
54: /* options for BoomerAMG */
55: PetscBool printstatistics;
57: /* options for BoomerAMG */
58: PetscInt cycletype;
59: PetscInt maxlevels;
60: PetscReal strongthreshold;
61: PetscReal maxrowsum;
62: PetscInt gridsweeps[3];
63: PetscInt coarsentype;
64: PetscInt measuretype;
65: PetscInt smoothtype;
66: PetscInt smoothnumlevels;
67: PetscInt eu_level; /* Number of levels for ILU(k) in Euclid */
68: PetscReal eu_droptolerance; /* Drop tolerance for ILU(k) in Euclid */
69: PetscInt eu_bj; /* Defines use of Block Jacobi ILU in Euclid */
70: PetscInt relaxtype[3];
71: PetscReal relaxweight;
72: PetscReal outerrelaxweight;
73: PetscInt relaxorder;
74: PetscReal truncfactor;
75: PetscBool applyrichardson;
76: PetscInt pmax;
77: PetscInt interptype;
78: PetscInt maxc;
79: PetscInt minc;
80: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
81: char *spgemm_type; // this is a global hypre parameter but is closely associated with BoomerAMG
82: #endif
83: /* GPU */
84: PetscBool keeptranspose;
85: PetscInt rap2;
86: PetscInt mod_rap2;
88: /* AIR */
89: PetscInt Rtype;
90: PetscReal Rstrongthreshold;
91: PetscReal Rfilterthreshold;
92: PetscInt Adroptype;
93: PetscReal Adroptol;
95: PetscInt agg_nl;
96: PetscInt agg_interptype;
97: PetscInt agg_num_paths;
98: PetscBool nodal_relax;
99: PetscInt nodal_relax_levels;
101: PetscInt nodal_coarsening;
102: PetscInt nodal_coarsening_diag;
103: PetscInt vec_interp_variant;
104: PetscInt vec_interp_qmax;
105: PetscBool vec_interp_smooth;
106: PetscInt interp_refine;
108: /* NearNullSpace support */
109: VecHYPRE_IJVector *hmnull;
110: HYPRE_ParVector *phmnull;
111: PetscInt n_hmnull;
112: Vec hmnull_constant;
114: /* options for AS (Auxiliary Space preconditioners) */
115: PetscInt as_print;
116: PetscInt as_max_iter;
117: PetscReal as_tol;
118: PetscInt as_relax_type;
119: PetscInt as_relax_times;
120: PetscReal as_relax_weight;
121: PetscReal as_omega;
122: PetscInt as_amg_alpha_opts[5]; /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for vector Poisson (AMS) or Curl problem (ADS) */
123: PetscReal as_amg_alpha_theta; /* AMG strength for vector Poisson (AMS) or Curl problem (ADS) */
124: PetscInt as_amg_beta_opts[5]; /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for scalar Poisson (AMS) or vector Poisson (ADS) */
125: PetscReal as_amg_beta_theta; /* AMG strength for scalar Poisson (AMS) or vector Poisson (ADS) */
126: PetscInt ams_cycle_type;
127: PetscInt ads_cycle_type;
129: /* additional data */
130: Mat G; /* MatHYPRE */
131: Mat C; /* MatHYPRE */
132: Mat alpha_Poisson; /* MatHYPRE */
133: Mat beta_Poisson; /* MatHYPRE */
135: /* extra information for AMS */
136: PetscInt dim; /* geometrical dimension */
137: VecHYPRE_IJVector coords[3];
138: VecHYPRE_IJVector constants[3];
139: VecHYPRE_IJVector interior;
140: Mat RT_PiFull, RT_Pi[3];
141: Mat ND_PiFull, ND_Pi[3];
142: PetscBool ams_beta_is_zero;
143: PetscBool ams_beta_is_zero_part;
144: PetscInt ams_proj_freq;
145: } PC_HYPRE;
147: PetscErrorCode PCHYPREGetSolver(PC pc, HYPRE_Solver *hsolver)
148: {
149: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
151: PetscFunctionBegin;
152: *hsolver = jac->hsolver;
153: PetscFunctionReturn(PETSC_SUCCESS);
154: }
156: /*
157: Matrices with AIJ format are created IN PLACE with using (I,J,data) from BoomerAMG. Since the data format in hypre_ParCSRMatrix
158: is different from that used in PETSc, the original hypre_ParCSRMatrix can not be used any more after call this routine.
159: It is used in PCHMG. Other users should avoid using this function.
160: */
161: static PetscErrorCode PCGetCoarseOperators_BoomerAMG(PC pc, PetscInt *nlevels, Mat *operators[])
162: {
163: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
164: PetscBool same = PETSC_FALSE;
165: PetscInt num_levels, l;
166: Mat *mattmp;
167: hypre_ParCSRMatrix **A_array;
169: PetscFunctionBegin;
170: PetscCall(PetscStrcmp(jac->hypre_type, "boomeramg", &same));
171: PetscCheck(same, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_NOTSAMETYPE, "Hypre type is not BoomerAMG ");
172: num_levels = hypre_ParAMGDataNumLevels((hypre_ParAMGData *)(jac->hsolver));
173: PetscCall(PetscMalloc1(num_levels, &mattmp));
174: A_array = hypre_ParAMGDataAArray((hypre_ParAMGData *)(jac->hsolver));
175: for (l = 1; l < num_levels; l++) {
176: PetscCall(MatCreateFromParCSR(A_array[l], MATAIJ, PETSC_OWN_POINTER, &(mattmp[num_levels - 1 - l])));
177: /* We want to own the data, and HYPRE can not touch this matrix any more */
178: A_array[l] = NULL;
179: }
180: *nlevels = num_levels;
181: *operators = mattmp;
182: PetscFunctionReturn(PETSC_SUCCESS);
183: }
185: /*
186: Matrices with AIJ format are created IN PLACE with using (I,J,data) from BoomerAMG. Since the data format in hypre_ParCSRMatrix
187: is different from that used in PETSc, the original hypre_ParCSRMatrix can not be used any more after call this routine.
188: It is used in PCHMG. Other users should avoid using this function.
189: */
190: static PetscErrorCode PCGetInterpolations_BoomerAMG(PC pc, PetscInt *nlevels, Mat *interpolations[])
191: {
192: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
193: PetscBool same = PETSC_FALSE;
194: PetscInt num_levels, l;
195: Mat *mattmp;
196: hypre_ParCSRMatrix **P_array;
198: PetscFunctionBegin;
199: PetscCall(PetscStrcmp(jac->hypre_type, "boomeramg", &same));
200: PetscCheck(same, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_NOTSAMETYPE, "Hypre type is not BoomerAMG ");
201: num_levels = hypre_ParAMGDataNumLevels((hypre_ParAMGData *)(jac->hsolver));
202: PetscCall(PetscMalloc1(num_levels, &mattmp));
203: P_array = hypre_ParAMGDataPArray((hypre_ParAMGData *)(jac->hsolver));
204: for (l = 1; l < num_levels; l++) {
205: PetscCall(MatCreateFromParCSR(P_array[num_levels - 1 - l], MATAIJ, PETSC_OWN_POINTER, &(mattmp[l - 1])));
206: /* We want to own the data, and HYPRE can not touch this matrix any more */
207: P_array[num_levels - 1 - l] = NULL;
208: }
209: *nlevels = num_levels;
210: *interpolations = mattmp;
211: PetscFunctionReturn(PETSC_SUCCESS);
212: }
214: /* Resets (frees) Hypre's representation of the near null space */
215: static PetscErrorCode PCHYPREResetNearNullSpace_Private(PC pc)
216: {
217: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
218: PetscInt i;
220: PetscFunctionBegin;
221: for (i = 0; i < jac->n_hmnull; i++) PetscCall(VecHYPRE_IJVectorDestroy(&jac->hmnull[i]));
222: PetscCall(PetscFree(jac->hmnull));
223: PetscCall(PetscFree(jac->phmnull));
224: PetscCall(VecDestroy(&jac->hmnull_constant));
225: jac->n_hmnull = 0;
226: PetscFunctionReturn(PETSC_SUCCESS);
227: }
229: static PetscErrorCode PCSetUp_HYPRE(PC pc)
230: {
231: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
232: Mat_HYPRE *hjac;
233: HYPRE_ParCSRMatrix hmat;
234: HYPRE_ParVector bv, xv;
235: PetscBool ishypre;
237: PetscFunctionBegin;
238: if (!jac->hypre_type) PetscCall(PCHYPRESetType(pc, "boomeramg"));
240: PetscCall(PetscObjectTypeCompare((PetscObject)pc->pmat, MATHYPRE, &ishypre));
241: if (!ishypre) {
242: PetscCall(MatDestroy(&jac->hpmat));
243: PetscCall(MatConvert(pc->pmat, MATHYPRE, MAT_INITIAL_MATRIX, &jac->hpmat));
244: } else {
245: PetscCall(PetscObjectReference((PetscObject)pc->pmat));
246: PetscCall(MatDestroy(&jac->hpmat));
247: jac->hpmat = pc->pmat;
248: }
249: /* allow debug */
250: PetscCall(MatViewFromOptions(jac->hpmat, NULL, "-pc_hypre_mat_view"));
251: hjac = (Mat_HYPRE *)(jac->hpmat->data);
253: /* special case for BoomerAMG */
254: if (jac->setup == HYPRE_BoomerAMGSetup) {
255: MatNullSpace mnull;
256: PetscBool has_const;
257: PetscInt bs, nvec, i;
258: const Vec *vecs;
260: PetscCall(MatGetBlockSize(pc->pmat, &bs));
261: if (bs > 1) PetscCallExternal(HYPRE_BoomerAMGSetNumFunctions, jac->hsolver, bs);
262: PetscCall(MatGetNearNullSpace(pc->mat, &mnull));
263: if (mnull) {
264: PetscCall(PCHYPREResetNearNullSpace_Private(pc));
265: PetscCall(MatNullSpaceGetVecs(mnull, &has_const, &nvec, &vecs));
266: PetscCall(PetscMalloc1(nvec + 1, &jac->hmnull));
267: PetscCall(PetscMalloc1(nvec + 1, &jac->phmnull));
268: for (i = 0; i < nvec; i++) {
269: PetscCall(VecHYPRE_IJVectorCreate(vecs[i]->map, &jac->hmnull[i]));
270: PetscCall(VecHYPRE_IJVectorCopy(vecs[i], jac->hmnull[i]));
271: PetscCallExternal(HYPRE_IJVectorGetObject, jac->hmnull[i]->ij, (void **)&jac->phmnull[i]);
272: }
273: if (has_const) {
274: PetscCall(MatCreateVecs(pc->pmat, &jac->hmnull_constant, NULL));
275: PetscCall(VecSet(jac->hmnull_constant, 1));
276: PetscCall(VecNormalize(jac->hmnull_constant, NULL));
277: PetscCall(VecHYPRE_IJVectorCreate(jac->hmnull_constant->map, &jac->hmnull[nvec]));
278: PetscCall(VecHYPRE_IJVectorCopy(jac->hmnull_constant, jac->hmnull[nvec]));
279: PetscCallExternal(HYPRE_IJVectorGetObject, jac->hmnull[nvec]->ij, (void **)&jac->phmnull[nvec]);
280: nvec++;
281: }
282: PetscCallExternal(HYPRE_BoomerAMGSetInterpVectors, jac->hsolver, nvec, jac->phmnull);
283: jac->n_hmnull = nvec;
284: }
285: }
287: /* special case for AMS */
288: if (jac->setup == HYPRE_AMSSetup) {
289: Mat_HYPRE *hm;
290: HYPRE_ParCSRMatrix parcsr;
291: if (!jac->coords[0] && !jac->constants[0] && !(jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1]))) {
292: 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()");
293: }
294: if (jac->dim) PetscCallExternal(HYPRE_AMSSetDimension, jac->hsolver, jac->dim);
295: if (jac->constants[0]) {
296: HYPRE_ParVector ozz, zoz, zzo = NULL;
297: PetscCallExternal(HYPRE_IJVectorGetObject, jac->constants[0]->ij, (void **)(&ozz));
298: PetscCallExternal(HYPRE_IJVectorGetObject, jac->constants[1]->ij, (void **)(&zoz));
299: if (jac->constants[2]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->constants[2]->ij, (void **)(&zzo));
300: PetscCallExternal(HYPRE_AMSSetEdgeConstantVectors, jac->hsolver, ozz, zoz, zzo);
301: }
302: if (jac->coords[0]) {
303: HYPRE_ParVector coords[3];
304: coords[0] = NULL;
305: coords[1] = NULL;
306: coords[2] = NULL;
307: if (jac->coords[0]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[0]->ij, (void **)(&coords[0]));
308: if (jac->coords[1]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[1]->ij, (void **)(&coords[1]));
309: if (jac->coords[2]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[2]->ij, (void **)(&coords[2]));
310: PetscCallExternal(HYPRE_AMSSetCoordinateVectors, jac->hsolver, coords[0], coords[1], coords[2]);
311: }
312: PetscCheck(jac->G, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "HYPRE AMS preconditioner needs the discrete gradient operator via PCHYPRESetDiscreteGradient");
313: hm = (Mat_HYPRE *)(jac->G->data);
314: PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&parcsr));
315: PetscCallExternal(HYPRE_AMSSetDiscreteGradient, jac->hsolver, parcsr);
316: if (jac->alpha_Poisson) {
317: hm = (Mat_HYPRE *)(jac->alpha_Poisson->data);
318: PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&parcsr));
319: PetscCallExternal(HYPRE_AMSSetAlphaPoissonMatrix, jac->hsolver, parcsr);
320: }
321: if (jac->ams_beta_is_zero) {
322: PetscCallExternal(HYPRE_AMSSetBetaPoissonMatrix, jac->hsolver, NULL);
323: } else if (jac->beta_Poisson) {
324: hm = (Mat_HYPRE *)(jac->beta_Poisson->data);
325: PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&parcsr));
326: PetscCallExternal(HYPRE_AMSSetBetaPoissonMatrix, jac->hsolver, parcsr);
327: } else if (jac->ams_beta_is_zero_part) {
328: if (jac->interior) {
329: HYPRE_ParVector interior = NULL;
330: PetscCallExternal(HYPRE_IJVectorGetObject, jac->interior->ij, (void **)(&interior));
331: PetscCallExternal(HYPRE_AMSSetInteriorNodes, jac->hsolver, interior);
332: } else {
333: jac->ams_beta_is_zero_part = PETSC_FALSE;
334: }
335: }
336: if (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1])) {
337: PetscInt i;
338: HYPRE_ParCSRMatrix nd_parcsrfull, nd_parcsr[3];
339: if (jac->ND_PiFull) {
340: hm = (Mat_HYPRE *)(jac->ND_PiFull->data);
341: PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&nd_parcsrfull));
342: } else {
343: nd_parcsrfull = NULL;
344: }
345: for (i = 0; i < 3; ++i) {
346: if (jac->ND_Pi[i]) {
347: hm = (Mat_HYPRE *)(jac->ND_Pi[i]->data);
348: PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&nd_parcsr[i]));
349: } else {
350: nd_parcsr[i] = NULL;
351: }
352: }
353: PetscCallExternal(HYPRE_AMSSetInterpolations, jac->hsolver, nd_parcsrfull, nd_parcsr[0], nd_parcsr[1], nd_parcsr[2]);
354: }
355: }
356: /* special case for ADS */
357: if (jac->setup == HYPRE_ADSSetup) {
358: Mat_HYPRE *hm;
359: HYPRE_ParCSRMatrix parcsr;
360: 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])))) {
361: SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "HYPRE ADS preconditioner needs either the coordinate vectors via PCSetCoordinates() or the interpolation matrices via PCHYPRESetInterpolations");
362: } else PetscCheck(jac->coords[1] && jac->coords[2], PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "HYPRE ADS preconditioner has been designed for three dimensional problems! For two dimensional problems, use HYPRE AMS instead");
363: PetscCheck(jac->G, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "HYPRE ADS preconditioner needs the discrete gradient operator via PCHYPRESetDiscreteGradient");
364: PetscCheck(jac->C, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "HYPRE ADS preconditioner needs the discrete curl operator via PCHYPRESetDiscreteGradient");
365: if (jac->coords[0]) {
366: HYPRE_ParVector coords[3];
367: coords[0] = NULL;
368: coords[1] = NULL;
369: coords[2] = NULL;
370: if (jac->coords[0]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[0]->ij, (void **)(&coords[0]));
371: if (jac->coords[1]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[1]->ij, (void **)(&coords[1]));
372: if (jac->coords[2]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[2]->ij, (void **)(&coords[2]));
373: PetscCallExternal(HYPRE_ADSSetCoordinateVectors, jac->hsolver, coords[0], coords[1], coords[2]);
374: }
375: hm = (Mat_HYPRE *)(jac->G->data);
376: PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&parcsr));
377: PetscCallExternal(HYPRE_ADSSetDiscreteGradient, jac->hsolver, parcsr);
378: hm = (Mat_HYPRE *)(jac->C->data);
379: PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&parcsr));
380: PetscCallExternal(HYPRE_ADSSetDiscreteCurl, jac->hsolver, parcsr);
381: if ((jac->RT_PiFull || (jac->RT_Pi[0] && jac->RT_Pi[1])) && (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1]))) {
382: PetscInt i;
383: HYPRE_ParCSRMatrix rt_parcsrfull, rt_parcsr[3];
384: HYPRE_ParCSRMatrix nd_parcsrfull, nd_parcsr[3];
385: if (jac->RT_PiFull) {
386: hm = (Mat_HYPRE *)(jac->RT_PiFull->data);
387: PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&rt_parcsrfull));
388: } else {
389: rt_parcsrfull = NULL;
390: }
391: for (i = 0; i < 3; ++i) {
392: if (jac->RT_Pi[i]) {
393: hm = (Mat_HYPRE *)(jac->RT_Pi[i]->data);
394: PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&rt_parcsr[i]));
395: } else {
396: rt_parcsr[i] = NULL;
397: }
398: }
399: if (jac->ND_PiFull) {
400: hm = (Mat_HYPRE *)(jac->ND_PiFull->data);
401: PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&nd_parcsrfull));
402: } else {
403: nd_parcsrfull = NULL;
404: }
405: for (i = 0; i < 3; ++i) {
406: if (jac->ND_Pi[i]) {
407: hm = (Mat_HYPRE *)(jac->ND_Pi[i]->data);
408: PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&nd_parcsr[i]));
409: } else {
410: nd_parcsr[i] = NULL;
411: }
412: }
413: PetscCallExternal(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]);
414: }
415: }
416: PetscCallExternal(HYPRE_IJMatrixGetObject, hjac->ij, (void **)&hmat);
417: PetscCallExternal(HYPRE_IJVectorGetObject, hjac->b->ij, (void **)&bv);
418: PetscCallExternal(HYPRE_IJVectorGetObject, hjac->x->ij, (void **)&xv);
419: PetscCall(PetscFPTrapPush(PETSC_FP_TRAP_OFF));
420: PetscCallExternal(jac->setup, jac->hsolver, hmat, bv, xv);
421: PetscCall(PetscFPTrapPop());
422: PetscFunctionReturn(PETSC_SUCCESS);
423: }
425: static PetscErrorCode PCApply_HYPRE(PC pc, Vec b, Vec x)
426: {
427: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
428: Mat_HYPRE *hjac = (Mat_HYPRE *)(jac->hpmat->data);
429: HYPRE_ParCSRMatrix hmat;
430: HYPRE_ParVector jbv, jxv;
432: PetscFunctionBegin;
433: PetscCall(PetscCitationsRegister(hypreCitation, &cite));
434: if (!jac->applyrichardson) PetscCall(VecSet(x, 0.0));
435: PetscCall(VecHYPRE_IJVectorPushVecRead(hjac->b, b));
436: if (jac->applyrichardson) PetscCall(VecHYPRE_IJVectorPushVec(hjac->x, x));
437: else PetscCall(VecHYPRE_IJVectorPushVecWrite(hjac->x, x));
438: PetscCallExternal(HYPRE_IJMatrixGetObject, hjac->ij, (void **)&hmat);
439: PetscCallExternal(HYPRE_IJVectorGetObject, hjac->b->ij, (void **)&jbv);
440: PetscCallExternal(HYPRE_IJVectorGetObject, hjac->x->ij, (void **)&jxv);
441: PetscStackCallExternalVoid(
442: "Hypre solve", do {
443: HYPRE_Int hierr = (*jac->solve)(jac->hsolver, hmat, jbv, jxv);
444: if (hierr) {
445: PetscCheck(hierr == HYPRE_ERROR_CONV, PETSC_COMM_SELF, PETSC_ERR_LIB, "Error in HYPRE solver, error code %d", (int)hierr);
446: hypre__global_error = 0;
447: }
448: } while (0));
450: if (jac->setup == HYPRE_AMSSetup && jac->ams_beta_is_zero_part) PetscCallExternal(HYPRE_AMSProjectOutGradients, jac->hsolver, jxv);
451: PetscCall(VecHYPRE_IJVectorPopVec(hjac->x));
452: PetscCall(VecHYPRE_IJVectorPopVec(hjac->b));
453: PetscFunctionReturn(PETSC_SUCCESS);
454: }
456: static PetscErrorCode PCReset_HYPRE(PC pc)
457: {
458: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
460: PetscFunctionBegin;
461: PetscCall(MatDestroy(&jac->hpmat));
462: PetscCall(MatDestroy(&jac->G));
463: PetscCall(MatDestroy(&jac->C));
464: PetscCall(MatDestroy(&jac->alpha_Poisson));
465: PetscCall(MatDestroy(&jac->beta_Poisson));
466: PetscCall(MatDestroy(&jac->RT_PiFull));
467: PetscCall(MatDestroy(&jac->RT_Pi[0]));
468: PetscCall(MatDestroy(&jac->RT_Pi[1]));
469: PetscCall(MatDestroy(&jac->RT_Pi[2]));
470: PetscCall(MatDestroy(&jac->ND_PiFull));
471: PetscCall(MatDestroy(&jac->ND_Pi[0]));
472: PetscCall(MatDestroy(&jac->ND_Pi[1]));
473: PetscCall(MatDestroy(&jac->ND_Pi[2]));
474: PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[0]));
475: PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[1]));
476: PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[2]));
477: PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[0]));
478: PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[1]));
479: PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[2]));
480: PetscCall(VecHYPRE_IJVectorDestroy(&jac->interior));
481: PetscCall(PCHYPREResetNearNullSpace_Private(pc));
482: jac->ams_beta_is_zero = PETSC_FALSE;
483: jac->ams_beta_is_zero_part = PETSC_FALSE;
484: jac->dim = 0;
485: PetscFunctionReturn(PETSC_SUCCESS);
486: }
488: static PetscErrorCode PCDestroy_HYPRE(PC pc)
489: {
490: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
492: PetscFunctionBegin;
493: PetscCall(PCReset_HYPRE(pc));
494: if (jac->destroy) PetscCallExternal(jac->destroy, jac->hsolver);
495: PetscCall(PetscFree(jac->hypre_type));
496: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
497: PetscCall(PetscFree(jac->spgemm_type));
498: #endif
499: if (jac->comm_hypre != MPI_COMM_NULL) PetscCall(PetscCommRestoreComm(PetscObjectComm((PetscObject)pc), &jac->comm_hypre));
500: PetscCall(PetscFree(pc->data));
502: PetscCall(PetscObjectChangeTypeName((PetscObject)pc, 0));
503: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetType_C", NULL));
504: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPREGetType_C", NULL));
505: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetDiscreteGradient_C", NULL));
506: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetDiscreteCurl_C", NULL));
507: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetInterpolations_C", NULL));
508: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetConstantEdgeVectors_C", NULL));
509: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetPoissonMatrix_C", NULL));
510: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetEdgeConstantVectors_C", NULL));
511: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPREAMSSetInteriorNodes_C", NULL));
512: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCGetInterpolations_C", NULL));
513: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCGetCoarseOperators_C", NULL));
514: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCMGGalerkinSetMatProductAlgorithm_C", NULL));
515: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCMGGalerkinGetMatProductAlgorithm_C", NULL));
516: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSetCoordinates_C", NULL));
517: PetscFunctionReturn(PETSC_SUCCESS);
518: }
520: static PetscErrorCode PCSetFromOptions_HYPRE_Pilut(PC pc, PetscOptionItems *PetscOptionsObject)
521: {
522: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
523: PetscBool flag;
525: PetscFunctionBegin;
526: PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE Pilut Options");
527: PetscCall(PetscOptionsInt("-pc_hypre_pilut_maxiter", "Number of iterations", "None", jac->maxiter, &jac->maxiter, &flag));
528: if (flag) PetscCallExternal(HYPRE_ParCSRPilutSetMaxIter, jac->hsolver, jac->maxiter);
529: PetscCall(PetscOptionsReal("-pc_hypre_pilut_tol", "Drop tolerance", "None", jac->tol, &jac->tol, &flag));
530: if (flag) PetscCallExternal(HYPRE_ParCSRPilutSetDropTolerance, jac->hsolver, jac->tol);
531: PetscCall(PetscOptionsInt("-pc_hypre_pilut_factorrowsize", "FactorRowSize", "None", jac->factorrowsize, &jac->factorrowsize, &flag));
532: if (flag) PetscCallExternal(HYPRE_ParCSRPilutSetFactorRowSize, jac->hsolver, jac->factorrowsize);
533: PetscOptionsHeadEnd();
534: PetscFunctionReturn(PETSC_SUCCESS);
535: }
537: static PetscErrorCode PCView_HYPRE_Pilut(PC pc, PetscViewer viewer)
538: {
539: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
540: PetscBool iascii;
542: PetscFunctionBegin;
543: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
544: if (iascii) {
545: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE Pilut preconditioning\n"));
546: if (jac->maxiter != PETSC_DEFAULT) {
547: PetscCall(PetscViewerASCIIPrintf(viewer, " maximum number of iterations %" PetscInt_FMT "\n", jac->maxiter));
548: } else {
549: PetscCall(PetscViewerASCIIPrintf(viewer, " default maximum number of iterations \n"));
550: }
551: if (jac->tol != PETSC_DEFAULT) {
552: PetscCall(PetscViewerASCIIPrintf(viewer, " drop tolerance %g\n", (double)jac->tol));
553: } else {
554: PetscCall(PetscViewerASCIIPrintf(viewer, " default drop tolerance \n"));
555: }
556: if (jac->factorrowsize != PETSC_DEFAULT) {
557: PetscCall(PetscViewerASCIIPrintf(viewer, " factor row size %" PetscInt_FMT "\n", jac->factorrowsize));
558: } else {
559: PetscCall(PetscViewerASCIIPrintf(viewer, " default factor row size \n"));
560: }
561: }
562: PetscFunctionReturn(PETSC_SUCCESS);
563: }
565: static PetscErrorCode PCSetFromOptions_HYPRE_Euclid(PC pc, PetscOptionItems *PetscOptionsObject)
566: {
567: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
568: PetscBool flag, eu_bj = jac->eu_bj ? PETSC_TRUE : PETSC_FALSE;
570: PetscFunctionBegin;
571: PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE Euclid Options");
572: PetscCall(PetscOptionsInt("-pc_hypre_euclid_level", "Factorization levels", "None", jac->eu_level, &jac->eu_level, &flag));
573: if (flag) PetscCallExternal(HYPRE_EuclidSetLevel, jac->hsolver, jac->eu_level);
575: PetscCall(PetscOptionsReal("-pc_hypre_euclid_droptolerance", "Drop tolerance for ILU(k) in Euclid", "None", jac->eu_droptolerance, &jac->eu_droptolerance, &flag));
576: if (flag) {
577: PetscMPIInt size;
579: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)pc), &size));
580: PetscCheck(size == 1, PetscObjectComm((PetscObject)pc), PETSC_ERR_SUP, "hypre's Euclid does not support a parallel drop tolerance");
581: PetscCallExternal(HYPRE_EuclidSetILUT, jac->hsolver, jac->eu_droptolerance);
582: }
584: PetscCall(PetscOptionsBool("-pc_hypre_euclid_bj", "Use Block Jacobi for ILU in Euclid", "None", eu_bj, &eu_bj, &flag));
585: if (flag) {
586: jac->eu_bj = eu_bj ? 1 : 0;
587: PetscCallExternal(HYPRE_EuclidSetBJ, jac->hsolver, jac->eu_bj);
588: }
589: PetscOptionsHeadEnd();
590: PetscFunctionReturn(PETSC_SUCCESS);
591: }
593: static PetscErrorCode PCView_HYPRE_Euclid(PC pc, PetscViewer viewer)
594: {
595: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
596: PetscBool iascii;
598: PetscFunctionBegin;
599: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
600: if (iascii) {
601: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE Euclid preconditioning\n"));
602: if (jac->eu_level != PETSC_DEFAULT) {
603: PetscCall(PetscViewerASCIIPrintf(viewer, " factorization levels %" PetscInt_FMT "\n", jac->eu_level));
604: } else {
605: PetscCall(PetscViewerASCIIPrintf(viewer, " default factorization levels \n"));
606: }
607: PetscCall(PetscViewerASCIIPrintf(viewer, " drop tolerance %g\n", (double)jac->eu_droptolerance));
608: PetscCall(PetscViewerASCIIPrintf(viewer, " use Block-Jacobi? %" PetscInt_FMT "\n", jac->eu_bj));
609: }
610: PetscFunctionReturn(PETSC_SUCCESS);
611: }
613: static PetscErrorCode PCApplyTranspose_HYPRE_BoomerAMG(PC pc, Vec b, Vec x)
614: {
615: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
616: Mat_HYPRE *hjac = (Mat_HYPRE *)(jac->hpmat->data);
617: HYPRE_ParCSRMatrix hmat;
618: HYPRE_ParVector jbv, jxv;
620: PetscFunctionBegin;
621: PetscCall(PetscCitationsRegister(hypreCitation, &cite));
622: PetscCall(VecSet(x, 0.0));
623: PetscCall(VecHYPRE_IJVectorPushVecRead(hjac->x, b));
624: PetscCall(VecHYPRE_IJVectorPushVecWrite(hjac->b, x));
626: PetscCallExternal(HYPRE_IJMatrixGetObject, hjac->ij, (void **)&hmat);
627: PetscCallExternal(HYPRE_IJVectorGetObject, hjac->b->ij, (void **)&jbv);
628: PetscCallExternal(HYPRE_IJVectorGetObject, hjac->x->ij, (void **)&jxv);
630: PetscStackCallExternalVoid(
631: "Hypre Transpose solve", do {
632: HYPRE_Int hierr = HYPRE_BoomerAMGSolveT(jac->hsolver, hmat, jbv, jxv);
633: if (hierr) {
634: /* error code of 1 in BoomerAMG merely means convergence not achieved */
635: PetscCheck(hierr == 1, PETSC_COMM_SELF, PETSC_ERR_LIB, "Error in HYPRE solver, error code %d", (int)hierr);
636: hypre__global_error = 0;
637: }
638: } while (0));
640: PetscCall(VecHYPRE_IJVectorPopVec(hjac->x));
641: PetscCall(VecHYPRE_IJVectorPopVec(hjac->b));
642: PetscFunctionReturn(PETSC_SUCCESS);
643: }
645: static PetscErrorCode PCMGGalerkinSetMatProductAlgorithm_HYPRE_BoomerAMG(PC pc, const char name[])
646: {
647: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
648: PetscBool flag;
650: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
651: PetscFunctionBegin;
652: if (jac->spgemm_type) {
653: PetscCall(PetscStrcmp(jac->spgemm_type, name, &flag));
654: PetscCheck(flag, PetscObjectComm((PetscObject)pc), PETSC_ERR_ORDER, "Cannot reset the HYPRE SpGEMM (really we can)");
655: PetscFunctionReturn(PETSC_SUCCESS);
656: } else {
657: PetscCall(PetscStrallocpy(name, &jac->spgemm_type));
658: }
659: PetscCall(PetscStrcmp("cusparse", jac->spgemm_type, &flag));
660: if (flag) {
661: PetscCallExternal(HYPRE_SetSpGemmUseCusparse, 1);
662: PetscFunctionReturn(PETSC_SUCCESS);
663: }
664: PetscCall(PetscStrcmp("hypre", jac->spgemm_type, &flag));
665: if (flag) {
666: PetscCallExternal(HYPRE_SetSpGemmUseCusparse, 0);
667: PetscFunctionReturn(PETSC_SUCCESS);
668: }
669: jac->spgemm_type = NULL;
670: SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown HYPRE SpGEM type %s; Choices are cusparse, hypre", name);
671: #endif
672: }
674: static PetscErrorCode PCMGGalerkinGetMatProductAlgorithm_HYPRE_BoomerAMG(PC pc, const char *spgemm[])
675: {
676: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
678: PetscFunctionBegin;
680: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
681: *spgemm = jac->spgemm_type;
682: #endif
683: PetscFunctionReturn(PETSC_SUCCESS);
684: }
686: static const char *HYPREBoomerAMGCycleType[] = {"", "V", "W"};
687: static const char *HYPREBoomerAMGCoarsenType[] = {"CLJP", "Ruge-Stueben", "", "modifiedRuge-Stueben", "", "", "Falgout", "", "PMIS", "", "HMIS"};
688: static const char *HYPREBoomerAMGMeasureType[] = {"local", "global"};
689: /* The following corresponds to HYPRE_BoomerAMGSetRelaxType which has many missing numbers in the enum */
690: static const char *HYPREBoomerAMGSmoothType[] = {"Schwarz-smoothers", "Pilut", "ParaSails", "Euclid"};
691: static const char *HYPREBoomerAMGRelaxType[] = {"Jacobi", "sequential-Gauss-Seidel", "seqboundary-Gauss-Seidel", "SOR/Jacobi", "backward-SOR/Jacobi", "" /* [5] hybrid chaotic Gauss-Seidel (works only with OpenMP) */, "symmetric-SOR/Jacobi", "" /* 7 */, "l1scaled-SOR/Jacobi", "Gaussian-elimination", "" /* 10 */, "" /* 11 */, "" /* 12 */, "l1-Gauss-Seidel" /* nonsymmetric */, "backward-l1-Gauss-Seidel" /* nonsymmetric */, "CG" /* non-stationary */, "Chebyshev", "FCF-Jacobi", "l1scaled-Jacobi"};
692: static const char *HYPREBoomerAMGInterpType[] = {"classical", "", "", "direct", "multipass", "multipass-wts", "ext+i", "ext+i-cc", "standard", "standard-wts", "block", "block-wtd", "FF", "FF1", "ext", "ad-wts", "ext-mm", "ext+i-mm", "ext+e-mm"};
693: static PetscErrorCode PCSetFromOptions_HYPRE_BoomerAMG(PC pc, PetscOptionItems *PetscOptionsObject)
694: {
695: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
696: PetscInt bs, n, indx, level;
697: PetscBool flg, tmp_truth;
698: double tmpdbl, twodbl[2];
699: const char *symtlist[] = {"nonsymmetric", "SPD", "nonsymmetric,SPD"};
700: const char *PCHYPRESpgemmTypes[] = {"cusparse", "hypre"};
702: PetscFunctionBegin;
703: PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE BoomerAMG Options");
704: PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_cycle_type", "Cycle type", "None", HYPREBoomerAMGCycleType + 1, 2, HYPREBoomerAMGCycleType[jac->cycletype], &indx, &flg));
705: if (flg) {
706: jac->cycletype = indx + 1;
707: PetscCallExternal(HYPRE_BoomerAMGSetCycleType, jac->hsolver, jac->cycletype);
708: }
709: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_max_levels", "Number of levels (of grids) allowed", "None", jac->maxlevels, &jac->maxlevels, &flg));
710: if (flg) {
711: PetscCheck(jac->maxlevels >= 2, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Number of levels %" PetscInt_FMT " must be at least two", jac->maxlevels);
712: PetscCallExternal(HYPRE_BoomerAMGSetMaxLevels, jac->hsolver, jac->maxlevels);
713: }
714: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_max_iter", "Maximum iterations used PER hypre call", "None", jac->maxiter, &jac->maxiter, &flg));
715: if (flg) {
716: PetscCheck(jac->maxiter >= 1, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Number of iterations %" PetscInt_FMT " must be at least one", jac->maxiter);
717: PetscCallExternal(HYPRE_BoomerAMGSetMaxIter, jac->hsolver, jac->maxiter);
718: }
719: PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_tol", "Convergence tolerance PER hypre call (0.0 = use a fixed number of iterations)", "None", jac->tol, &jac->tol, &flg));
720: if (flg) {
721: PetscCheck(jac->tol >= 0.0, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Tolerance %g must be greater than or equal to zero", (double)jac->tol);
722: PetscCallExternal(HYPRE_BoomerAMGSetTol, jac->hsolver, jac->tol);
723: }
724: bs = 1;
725: if (pc->pmat) PetscCall(MatGetBlockSize(pc->pmat, &bs));
726: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_numfunctions", "Number of functions", "HYPRE_BoomerAMGSetNumFunctions", bs, &bs, &flg));
727: if (flg) PetscCallExternal(HYPRE_BoomerAMGSetNumFunctions, jac->hsolver, bs);
729: PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_truncfactor", "Truncation factor for interpolation (0=no truncation)", "None", jac->truncfactor, &jac->truncfactor, &flg));
730: if (flg) {
731: PetscCheck(jac->truncfactor >= 0.0, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Truncation factor %g must be great than or equal zero", (double)jac->truncfactor);
732: PetscCallExternal(HYPRE_BoomerAMGSetTruncFactor, jac->hsolver, jac->truncfactor);
733: }
735: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_P_max", "Max elements per row for interpolation operator (0=unlimited)", "None", jac->pmax, &jac->pmax, &flg));
736: if (flg) {
737: PetscCheck(jac->pmax >= 0, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "P_max %" PetscInt_FMT " must be greater than or equal to zero", jac->pmax);
738: PetscCallExternal(HYPRE_BoomerAMGSetPMaxElmts, jac->hsolver, jac->pmax);
739: }
741: PetscCall(PetscOptionsRangeInt("-pc_hypre_boomeramg_agg_nl", "Number of levels of aggressive coarsening", "None", jac->agg_nl, &jac->agg_nl, &flg, 0, jac->maxlevels));
742: if (flg) PetscCallExternal(HYPRE_BoomerAMGSetAggNumLevels, jac->hsolver, jac->agg_nl);
744: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_agg_num_paths", "Number of paths for aggressive coarsening", "None", jac->agg_num_paths, &jac->agg_num_paths, &flg));
745: if (flg) {
746: PetscCheck(jac->agg_num_paths >= 1, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Number of paths %" PetscInt_FMT " must be greater than or equal to 1", jac->agg_num_paths);
747: PetscCallExternal(HYPRE_BoomerAMGSetNumPaths, jac->hsolver, jac->agg_num_paths);
748: }
750: PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_strong_threshold", "Threshold for being strongly connected", "None", jac->strongthreshold, &jac->strongthreshold, &flg));
751: if (flg) {
752: PetscCheck(jac->strongthreshold >= 0.0, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Strong threshold %g must be great than or equal zero", (double)jac->strongthreshold);
753: PetscCallExternal(HYPRE_BoomerAMGSetStrongThreshold, jac->hsolver, jac->strongthreshold);
754: }
755: PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_max_row_sum", "Maximum row sum", "None", jac->maxrowsum, &jac->maxrowsum, &flg));
756: if (flg) {
757: PetscCheck(jac->maxrowsum >= 0.0, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Maximum row sum %g must be greater than zero", (double)jac->maxrowsum);
758: PetscCheck(jac->maxrowsum <= 1.0, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Maximum row sum %g must be less than or equal one", (double)jac->maxrowsum);
759: PetscCallExternal(HYPRE_BoomerAMGSetMaxRowSum, jac->hsolver, jac->maxrowsum);
760: }
762: /* Grid sweeps */
763: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_all", "Number of sweeps for the up and down grid levels", "None", jac->gridsweeps[0], &indx, &flg));
764: if (flg) {
765: PetscCallExternal(HYPRE_BoomerAMGSetNumSweeps, jac->hsolver, indx);
766: /* modify the jac structure so we can view the updated options with PC_View */
767: jac->gridsweeps[0] = indx;
768: jac->gridsweeps[1] = indx;
769: /*defaults coarse to 1 */
770: jac->gridsweeps[2] = 1;
771: }
772: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_nodal_coarsen", "Use a nodal based coarsening 1-6", "HYPRE_BoomerAMGSetNodal", jac->nodal_coarsening, &jac->nodal_coarsening, &flg));
773: if (flg) PetscCallExternal(HYPRE_BoomerAMGSetNodal, jac->hsolver, jac->nodal_coarsening);
774: PetscCall(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));
775: if (flg) PetscCallExternal(HYPRE_BoomerAMGSetNodalDiag, jac->hsolver, jac->nodal_coarsening_diag);
776: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_vec_interp_variant", "Variant of algorithm 1-3", "HYPRE_BoomerAMGSetInterpVecVariant", jac->vec_interp_variant, &jac->vec_interp_variant, &flg));
777: if (flg) PetscCallExternal(HYPRE_BoomerAMGSetInterpVecVariant, jac->hsolver, jac->vec_interp_variant);
778: PetscCall(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));
779: if (flg) PetscCallExternal(HYPRE_BoomerAMGSetInterpVecQMax, jac->hsolver, jac->vec_interp_qmax);
780: PetscCall(PetscOptionsBool("-pc_hypre_boomeramg_vec_interp_smooth", "Whether to smooth the interpolation vectors", "HYPRE_BoomerAMGSetSmoothInterpVectors", jac->vec_interp_smooth, &jac->vec_interp_smooth, &flg));
781: if (flg) PetscCallExternal(HYPRE_BoomerAMGSetSmoothInterpVectors, jac->hsolver, jac->vec_interp_smooth);
782: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_interp_refine", "Preprocess the interpolation matrix through iterative weight refinement", "HYPRE_BoomerAMGSetInterpRefine", jac->interp_refine, &jac->interp_refine, &flg));
783: if (flg) PetscCallExternal(HYPRE_BoomerAMGSetInterpRefine, jac->hsolver, jac->interp_refine);
784: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_down", "Number of sweeps for the down cycles", "None", jac->gridsweeps[0], &indx, &flg));
785: if (flg) {
786: PetscCallExternal(HYPRE_BoomerAMGSetCycleNumSweeps, jac->hsolver, indx, 1);
787: jac->gridsweeps[0] = indx;
788: }
789: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_up", "Number of sweeps for the up cycles", "None", jac->gridsweeps[1], &indx, &flg));
790: if (flg) {
791: PetscCallExternal(HYPRE_BoomerAMGSetCycleNumSweeps, jac->hsolver, indx, 2);
792: jac->gridsweeps[1] = indx;
793: }
794: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_coarse", "Number of sweeps for the coarse level", "None", jac->gridsweeps[2], &indx, &flg));
795: if (flg) {
796: PetscCallExternal(HYPRE_BoomerAMGSetCycleNumSweeps, jac->hsolver, indx, 3);
797: jac->gridsweeps[2] = indx;
798: }
800: /* Smooth type */
801: PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_smooth_type", "Enable more complex smoothers", "None", HYPREBoomerAMGSmoothType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGSmoothType), HYPREBoomerAMGSmoothType[0], &indx, &flg));
802: if (flg) {
803: jac->smoothtype = indx;
804: PetscCallExternal(HYPRE_BoomerAMGSetSmoothType, jac->hsolver, indx + 6);
805: jac->smoothnumlevels = 25;
806: PetscCallExternal(HYPRE_BoomerAMGSetSmoothNumLevels, jac->hsolver, 25);
807: }
809: /* Number of smoothing levels */
810: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_smooth_num_levels", "Number of levels on which more complex smoothers are used", "None", 25, &indx, &flg));
811: if (flg && (jac->smoothtype != -1)) {
812: jac->smoothnumlevels = indx;
813: PetscCallExternal(HYPRE_BoomerAMGSetSmoothNumLevels, jac->hsolver, indx);
814: }
816: /* Number of levels for ILU(k) for Euclid */
817: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_eu_level", "Number of levels for ILU(k) in Euclid smoother", "None", 0, &indx, &flg));
818: if (flg && (jac->smoothtype == 3)) {
819: jac->eu_level = indx;
820: PetscCallExternal(HYPRE_BoomerAMGSetEuLevel, jac->hsolver, indx);
821: }
823: /* Filter for ILU(k) for Euclid */
824: double droptolerance;
825: PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_eu_droptolerance", "Drop tolerance for ILU(k) in Euclid smoother", "None", 0, &droptolerance, &flg));
826: if (flg && (jac->smoothtype == 3)) {
827: jac->eu_droptolerance = droptolerance;
828: PetscCallExternal(HYPRE_BoomerAMGSetEuLevel, jac->hsolver, droptolerance);
829: }
831: /* Use Block Jacobi ILUT for Euclid */
832: PetscCall(PetscOptionsBool("-pc_hypre_boomeramg_eu_bj", "Use Block Jacobi for ILU in Euclid smoother?", "None", PETSC_FALSE, &tmp_truth, &flg));
833: if (flg && (jac->smoothtype == 3)) {
834: jac->eu_bj = tmp_truth;
835: PetscCallExternal(HYPRE_BoomerAMGSetEuBJ, jac->hsolver, jac->eu_bj);
836: }
838: /* Relax type */
839: PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_relax_type_all", "Relax type for the up and down cycles", "None", HYPREBoomerAMGRelaxType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGRelaxType), HYPREBoomerAMGRelaxType[6], &indx, &flg));
840: if (flg) {
841: jac->relaxtype[0] = jac->relaxtype[1] = indx;
842: PetscCallExternal(HYPRE_BoomerAMGSetRelaxType, jac->hsolver, indx);
843: /* by default, coarse type set to 9 */
844: jac->relaxtype[2] = 9;
845: PetscCallExternal(HYPRE_BoomerAMGSetCycleRelaxType, jac->hsolver, 9, 3);
846: }
847: PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_relax_type_down", "Relax type for the down cycles", "None", HYPREBoomerAMGRelaxType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGRelaxType), HYPREBoomerAMGRelaxType[6], &indx, &flg));
848: if (flg) {
849: jac->relaxtype[0] = indx;
850: PetscCallExternal(HYPRE_BoomerAMGSetCycleRelaxType, jac->hsolver, indx, 1);
851: }
852: PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_relax_type_up", "Relax type for the up cycles", "None", HYPREBoomerAMGRelaxType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGRelaxType), HYPREBoomerAMGRelaxType[6], &indx, &flg));
853: if (flg) {
854: jac->relaxtype[1] = indx;
855: PetscCallExternal(HYPRE_BoomerAMGSetCycleRelaxType, jac->hsolver, indx, 2);
856: }
857: PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_relax_type_coarse", "Relax type on coarse grid", "None", HYPREBoomerAMGRelaxType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGRelaxType), HYPREBoomerAMGRelaxType[9], &indx, &flg));
858: if (flg) {
859: jac->relaxtype[2] = indx;
860: PetscCallExternal(HYPRE_BoomerAMGSetCycleRelaxType, jac->hsolver, indx, 3);
861: }
863: /* Relaxation Weight */
864: PetscCall(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));
865: if (flg) {
866: PetscCallExternal(HYPRE_BoomerAMGSetRelaxWt, jac->hsolver, tmpdbl);
867: jac->relaxweight = tmpdbl;
868: }
870: n = 2;
871: twodbl[0] = twodbl[1] = 1.0;
872: PetscCall(PetscOptionsRealArray("-pc_hypre_boomeramg_relax_weight_level", "Set the relaxation weight for a particular level (weight,level)", "None", twodbl, &n, &flg));
873: if (flg) {
874: PetscCheck(n == 2, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Relax weight level: you must provide 2 values separated by a comma (and no space), you provided %" PetscInt_FMT, n);
875: indx = (int)PetscAbsReal(twodbl[1]);
876: PetscCallExternal(HYPRE_BoomerAMGSetLevelRelaxWt, jac->hsolver, twodbl[0], indx);
877: }
879: /* Outer relaxation Weight */
880: PetscCall(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));
881: if (flg) {
882: PetscCallExternal(HYPRE_BoomerAMGSetOuterWt, jac->hsolver, tmpdbl);
883: jac->outerrelaxweight = tmpdbl;
884: }
886: n = 2;
887: twodbl[0] = twodbl[1] = 1.0;
888: PetscCall(PetscOptionsRealArray("-pc_hypre_boomeramg_outer_relax_weight_level", "Set the outer relaxation weight for a particular level (weight,level)", "None", twodbl, &n, &flg));
889: if (flg) {
890: PetscCheck(n == 2, 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 %" PetscInt_FMT, n);
891: indx = (int)PetscAbsReal(twodbl[1]);
892: PetscCallExternal(HYPRE_BoomerAMGSetLevelOuterWt, jac->hsolver, twodbl[0], indx);
893: }
895: /* the Relax Order */
896: PetscCall(PetscOptionsBool("-pc_hypre_boomeramg_no_CF", "Do not use CF-relaxation", "None", PETSC_FALSE, &tmp_truth, &flg));
898: if (flg && tmp_truth) {
899: jac->relaxorder = 0;
900: PetscCallExternal(HYPRE_BoomerAMGSetRelaxOrder, jac->hsolver, jac->relaxorder);
901: }
902: PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_measure_type", "Measure type", "None", HYPREBoomerAMGMeasureType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGMeasureType), HYPREBoomerAMGMeasureType[0], &indx, &flg));
903: if (flg) {
904: jac->measuretype = indx;
905: PetscCallExternal(HYPRE_BoomerAMGSetMeasureType, jac->hsolver, jac->measuretype);
906: }
907: /* update list length 3/07 */
908: PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_coarsen_type", "Coarsen type", "None", HYPREBoomerAMGCoarsenType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGCoarsenType), HYPREBoomerAMGCoarsenType[6], &indx, &flg));
909: if (flg) {
910: jac->coarsentype = indx;
911: PetscCallExternal(HYPRE_BoomerAMGSetCoarsenType, jac->hsolver, jac->coarsentype);
912: }
914: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_max_coarse_size", "Maximum size of coarsest grid", "None", jac->maxc, &jac->maxc, &flg));
915: if (flg) PetscCallExternal(HYPRE_BoomerAMGSetMaxCoarseSize, jac->hsolver, jac->maxc);
916: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_min_coarse_size", "Minimum size of coarsest grid", "None", jac->minc, &jac->minc, &flg));
917: if (flg) PetscCallExternal(HYPRE_BoomerAMGSetMinCoarseSize, jac->hsolver, jac->minc);
918: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
919: // global parameter but is closely associated with BoomerAMG
920: PetscCall(PetscOptionsEList("-pc_mg_galerkin_mat_product_algorithm", "Type of SpGEMM to use in hypre (only for now)", "PCMGGalerkinSetMatProductAlgorithm", PCHYPRESpgemmTypes, PETSC_STATIC_ARRAY_LENGTH(PCHYPRESpgemmTypes), PCHYPRESpgemmTypes[0], &indx, &flg));
921: if (!flg) indx = 0;
922: PetscCall(PCMGGalerkinSetMatProductAlgorithm_HYPRE_BoomerAMG(pc, PCHYPRESpgemmTypes[indx]));
923: #endif
924: /* AIR */
925: #if PETSC_PKG_HYPRE_VERSION_GE(2, 18, 0)
926: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_restriction_type", "Type of AIR method (distance 1 or 2, 0 means no AIR)", "None", jac->Rtype, &jac->Rtype, NULL));
927: PetscCallExternal(HYPRE_BoomerAMGSetRestriction, jac->hsolver, jac->Rtype);
928: if (jac->Rtype) {
929: jac->interptype = 100; /* no way we can pass this with strings... Set it as default as in MFEM, then users can still customize it back to a different one */
931: PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_strongthresholdR", "Threshold for R", "None", jac->Rstrongthreshold, &jac->Rstrongthreshold, NULL));
932: PetscCallExternal(HYPRE_BoomerAMGSetStrongThresholdR, jac->hsolver, jac->Rstrongthreshold);
934: PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_filterthresholdR", "Filter threshold for R", "None", jac->Rfilterthreshold, &jac->Rfilterthreshold, NULL));
935: PetscCallExternal(HYPRE_BoomerAMGSetFilterThresholdR, jac->hsolver, jac->Rfilterthreshold);
937: PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_Adroptol", "Defines the drop tolerance for the A-matrices from the 2nd level of AMG", "None", jac->Adroptol, &jac->Adroptol, NULL));
938: PetscCallExternal(HYPRE_BoomerAMGSetADropTol, jac->hsolver, jac->Adroptol);
940: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_Adroptype", "Drops the entries that are not on the diagonal and smaller than its row norm: type 1: 1-norm, 2: 2-norm, -1: infinity norm", "None", jac->Adroptype, &jac->Adroptype, NULL));
941: PetscCallExternal(HYPRE_BoomerAMGSetADropType, jac->hsolver, jac->Adroptype);
942: }
943: #endif
945: #if PETSC_PKG_HYPRE_VERSION_LE(9, 9, 9)
946: PetscCheck(!jac->Rtype || !jac->agg_nl, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_INCOMP, "-pc_hypre_boomeramg_restriction_type (%" PetscInt_FMT ") and -pc_hypre_boomeramg_agg_nl (%" PetscInt_FMT ")", jac->Rtype, jac->agg_nl);
947: #endif
949: /* new 3/07 */
950: PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_interp_type", "Interpolation type", "None", HYPREBoomerAMGInterpType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGInterpType), HYPREBoomerAMGInterpType[0], &indx, &flg));
951: if (flg || jac->Rtype) {
952: if (flg) jac->interptype = indx;
953: PetscCallExternal(HYPRE_BoomerAMGSetInterpType, jac->hsolver, jac->interptype);
954: }
956: PetscCall(PetscOptionsName("-pc_hypre_boomeramg_print_statistics", "Print statistics", "None", &flg));
957: if (flg) {
958: level = 3;
959: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_print_statistics", "Print statistics", "None", level, &level, NULL));
961: jac->printstatistics = PETSC_TRUE;
962: PetscCallExternal(HYPRE_BoomerAMGSetPrintLevel, jac->hsolver, level);
963: }
965: PetscCall(PetscOptionsName("-pc_hypre_boomeramg_print_debug", "Print debug information", "None", &flg));
966: if (flg) {
967: level = 3;
968: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_print_debug", "Print debug information", "None", level, &level, NULL));
970: jac->printstatistics = PETSC_TRUE;
971: PetscCallExternal(HYPRE_BoomerAMGSetDebugFlag, jac->hsolver, level);
972: }
974: PetscCall(PetscOptionsBool("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None", PETSC_FALSE, &tmp_truth, &flg));
975: if (flg && tmp_truth) {
976: PetscInt tmp_int;
977: PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None", jac->nodal_relax_levels, &tmp_int, &flg));
978: if (flg) jac->nodal_relax_levels = tmp_int;
979: PetscCallExternal(HYPRE_BoomerAMGSetSmoothType, jac->hsolver, 6);
980: PetscCallExternal(HYPRE_BoomerAMGSetDomainType, jac->hsolver, 1);
981: PetscCallExternal(HYPRE_BoomerAMGSetOverlap, jac->hsolver, 0);
982: PetscCallExternal(HYPRE_BoomerAMGSetSmoothNumLevels, jac->hsolver, jac->nodal_relax_levels);
983: }
985: PetscCall(PetscOptionsBool("-pc_hypre_boomeramg_keeptranspose", "Avoid transpose matvecs in preconditioner application", "None", jac->keeptranspose, &jac->keeptranspose, NULL));
986: PetscCallExternal(HYPRE_BoomerAMGSetKeepTranspose, jac->hsolver, jac->keeptranspose ? 1 : 0);
988: /* options for ParaSails solvers */
989: PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_parasails_sym", "Symmetry of matrix and preconditioner", "None", symtlist, PETSC_STATIC_ARRAY_LENGTH(symtlist), symtlist[0], &indx, &flg));
990: if (flg) {
991: jac->symt = indx;
992: PetscCallExternal(HYPRE_BoomerAMGSetSym, jac->hsolver, jac->symt);
993: }
995: PetscOptionsHeadEnd();
996: PetscFunctionReturn(PETSC_SUCCESS);
997: }
999: 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)
1000: {
1001: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1002: HYPRE_Int oits;
1004: PetscFunctionBegin;
1005: PetscCall(PetscCitationsRegister(hypreCitation, &cite));
1006: PetscCallExternal(HYPRE_BoomerAMGSetMaxIter, jac->hsolver, its * jac->maxiter);
1007: PetscCallExternal(HYPRE_BoomerAMGSetTol, jac->hsolver, rtol);
1008: jac->applyrichardson = PETSC_TRUE;
1009: PetscCall(PCApply_HYPRE(pc, b, y));
1010: jac->applyrichardson = PETSC_FALSE;
1011: PetscCallExternal(HYPRE_BoomerAMGGetNumIterations, jac->hsolver, &oits);
1012: *outits = oits;
1013: if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
1014: else *reason = PCRICHARDSON_CONVERGED_RTOL;
1015: PetscCallExternal(HYPRE_BoomerAMGSetTol, jac->hsolver, jac->tol);
1016: PetscCallExternal(HYPRE_BoomerAMGSetMaxIter, jac->hsolver, jac->maxiter);
1017: PetscFunctionReturn(PETSC_SUCCESS);
1018: }
1020: static PetscErrorCode PCView_HYPRE_BoomerAMG(PC pc, PetscViewer viewer)
1021: {
1022: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1023: PetscBool iascii;
1025: PetscFunctionBegin;
1026: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
1027: if (iascii) {
1028: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE BoomerAMG preconditioning\n"));
1029: PetscCall(PetscViewerASCIIPrintf(viewer, " Cycle type %s\n", HYPREBoomerAMGCycleType[jac->cycletype]));
1030: PetscCall(PetscViewerASCIIPrintf(viewer, " Maximum number of levels %" PetscInt_FMT "\n", jac->maxlevels));
1031: PetscCall(PetscViewerASCIIPrintf(viewer, " Maximum number of iterations PER hypre call %" PetscInt_FMT "\n", jac->maxiter));
1032: PetscCall(PetscViewerASCIIPrintf(viewer, " Convergence tolerance PER hypre call %g\n", (double)jac->tol));
1033: PetscCall(PetscViewerASCIIPrintf(viewer, " Threshold for strong coupling %g\n", (double)jac->strongthreshold));
1034: PetscCall(PetscViewerASCIIPrintf(viewer, " Interpolation truncation factor %g\n", (double)jac->truncfactor));
1035: PetscCall(PetscViewerASCIIPrintf(viewer, " Interpolation: max elements per row %" PetscInt_FMT "\n", jac->pmax));
1036: if (jac->interp_refine) PetscCall(PetscViewerASCIIPrintf(viewer, " Interpolation: number of steps of weighted refinement %" PetscInt_FMT "\n", jac->interp_refine));
1037: PetscCall(PetscViewerASCIIPrintf(viewer, " Number of levels of aggressive coarsening %" PetscInt_FMT "\n", jac->agg_nl));
1038: PetscCall(PetscViewerASCIIPrintf(viewer, " Number of paths for aggressive coarsening %" PetscInt_FMT "\n", jac->agg_num_paths));
1040: PetscCall(PetscViewerASCIIPrintf(viewer, " Maximum row sums %g\n", (double)jac->maxrowsum));
1042: PetscCall(PetscViewerASCIIPrintf(viewer, " Sweeps down %" PetscInt_FMT "\n", jac->gridsweeps[0]));
1043: PetscCall(PetscViewerASCIIPrintf(viewer, " Sweeps up %" PetscInt_FMT "\n", jac->gridsweeps[1]));
1044: PetscCall(PetscViewerASCIIPrintf(viewer, " Sweeps on coarse %" PetscInt_FMT "\n", jac->gridsweeps[2]));
1046: PetscCall(PetscViewerASCIIPrintf(viewer, " Relax down %s\n", HYPREBoomerAMGRelaxType[jac->relaxtype[0]]));
1047: PetscCall(PetscViewerASCIIPrintf(viewer, " Relax up %s\n", HYPREBoomerAMGRelaxType[jac->relaxtype[1]]));
1048: PetscCall(PetscViewerASCIIPrintf(viewer, " Relax on coarse %s\n", HYPREBoomerAMGRelaxType[jac->relaxtype[2]]));
1050: PetscCall(PetscViewerASCIIPrintf(viewer, " Relax weight (all) %g\n", (double)jac->relaxweight));
1051: PetscCall(PetscViewerASCIIPrintf(viewer, " Outer relax weight (all) %g\n", (double)jac->outerrelaxweight));
1053: if (jac->relaxorder) {
1054: PetscCall(PetscViewerASCIIPrintf(viewer, " Using CF-relaxation\n"));
1055: } else {
1056: PetscCall(PetscViewerASCIIPrintf(viewer, " Not using CF-relaxation\n"));
1057: }
1058: if (jac->smoothtype != -1) {
1059: PetscCall(PetscViewerASCIIPrintf(viewer, " Smooth type %s\n", HYPREBoomerAMGSmoothType[jac->smoothtype]));
1060: PetscCall(PetscViewerASCIIPrintf(viewer, " Smooth num levels %" PetscInt_FMT "\n", jac->smoothnumlevels));
1061: } else {
1062: PetscCall(PetscViewerASCIIPrintf(viewer, " Not using more complex smoothers.\n"));
1063: }
1064: if (jac->smoothtype == 3) {
1065: PetscCall(PetscViewerASCIIPrintf(viewer, " Euclid ILU(k) levels %" PetscInt_FMT "\n", jac->eu_level));
1066: PetscCall(PetscViewerASCIIPrintf(viewer, " Euclid ILU(k) drop tolerance %g\n", (double)jac->eu_droptolerance));
1067: PetscCall(PetscViewerASCIIPrintf(viewer, " Euclid ILU use Block-Jacobi? %" PetscInt_FMT "\n", jac->eu_bj));
1068: }
1069: PetscCall(PetscViewerASCIIPrintf(viewer, " Measure type %s\n", HYPREBoomerAMGMeasureType[jac->measuretype]));
1070: PetscCall(PetscViewerASCIIPrintf(viewer, " Coarsen type %s\n", HYPREBoomerAMGCoarsenType[jac->coarsentype]));
1071: PetscCall(PetscViewerASCIIPrintf(viewer, " Interpolation type %s\n", jac->interptype != 100 ? HYPREBoomerAMGInterpType[jac->interptype] : "1pt"));
1072: if (jac->nodal_coarsening) PetscCall(PetscViewerASCIIPrintf(viewer, " Using nodal coarsening with HYPRE_BOOMERAMGSetNodal() %" PetscInt_FMT "\n", jac->nodal_coarsening));
1073: if (jac->vec_interp_variant) {
1074: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE_BoomerAMGSetInterpVecVariant() %" PetscInt_FMT "\n", jac->vec_interp_variant));
1075: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE_BoomerAMGSetInterpVecQMax() %" PetscInt_FMT "\n", jac->vec_interp_qmax));
1076: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE_BoomerAMGSetSmoothInterpVectors() %d\n", jac->vec_interp_smooth));
1077: }
1078: if (jac->nodal_relax) PetscCall(PetscViewerASCIIPrintf(viewer, " Using nodal relaxation via Schwarz smoothing on levels %" PetscInt_FMT "\n", jac->nodal_relax_levels));
1079: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
1080: PetscCall(PetscViewerASCIIPrintf(viewer, " SpGEMM type %s\n", jac->spgemm_type));
1081: #endif
1082: /* AIR */
1083: if (jac->Rtype) {
1084: PetscCall(PetscViewerASCIIPrintf(viewer, " Using approximate ideal restriction type %" PetscInt_FMT "\n", jac->Rtype));
1085: PetscCall(PetscViewerASCIIPrintf(viewer, " Threshold for R %g\n", (double)jac->Rstrongthreshold));
1086: PetscCall(PetscViewerASCIIPrintf(viewer, " Filter for R %g\n", (double)jac->Rfilterthreshold));
1087: PetscCall(PetscViewerASCIIPrintf(viewer, " A drop tolerance %g\n", (double)jac->Adroptol));
1088: PetscCall(PetscViewerASCIIPrintf(viewer, " A drop type %" PetscInt_FMT "\n", jac->Adroptype));
1089: }
1090: }
1091: PetscFunctionReturn(PETSC_SUCCESS);
1092: }
1094: static PetscErrorCode PCSetFromOptions_HYPRE_ParaSails(PC pc, PetscOptionItems *PetscOptionsObject)
1095: {
1096: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1097: PetscInt indx;
1098: PetscBool flag;
1099: const char *symtlist[] = {"nonsymmetric", "SPD", "nonsymmetric,SPD"};
1101: PetscFunctionBegin;
1102: PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE ParaSails Options");
1103: PetscCall(PetscOptionsInt("-pc_hypre_parasails_nlevels", "Number of number of levels", "None", jac->nlevels, &jac->nlevels, 0));
1104: PetscCall(PetscOptionsReal("-pc_hypre_parasails_thresh", "Threshold", "None", jac->threshold, &jac->threshold, &flag));
1105: if (flag) PetscCallExternal(HYPRE_ParaSailsSetParams, jac->hsolver, jac->threshold, jac->nlevels);
1107: PetscCall(PetscOptionsReal("-pc_hypre_parasails_filter", "filter", "None", jac->filter, &jac->filter, &flag));
1108: if (flag) PetscCallExternal(HYPRE_ParaSailsSetFilter, jac->hsolver, jac->filter);
1110: PetscCall(PetscOptionsReal("-pc_hypre_parasails_loadbal", "Load balance", "None", jac->loadbal, &jac->loadbal, &flag));
1111: if (flag) PetscCallExternal(HYPRE_ParaSailsSetLoadbal, jac->hsolver, jac->loadbal);
1113: PetscCall(PetscOptionsBool("-pc_hypre_parasails_logging", "Print info to screen", "None", (PetscBool)jac->logging, (PetscBool *)&jac->logging, &flag));
1114: if (flag) PetscCallExternal(HYPRE_ParaSailsSetLogging, jac->hsolver, jac->logging);
1116: PetscCall(PetscOptionsBool("-pc_hypre_parasails_reuse", "Reuse nonzero pattern in preconditioner", "None", (PetscBool)jac->ruse, (PetscBool *)&jac->ruse, &flag));
1117: if (flag) PetscCallExternal(HYPRE_ParaSailsSetReuse, jac->hsolver, jac->ruse);
1119: PetscCall(PetscOptionsEList("-pc_hypre_parasails_sym", "Symmetry of matrix and preconditioner", "None", symtlist, PETSC_STATIC_ARRAY_LENGTH(symtlist), symtlist[0], &indx, &flag));
1120: if (flag) {
1121: jac->symt = indx;
1122: PetscCallExternal(HYPRE_ParaSailsSetSym, jac->hsolver, jac->symt);
1123: }
1125: PetscOptionsHeadEnd();
1126: PetscFunctionReturn(PETSC_SUCCESS);
1127: }
1129: static PetscErrorCode PCView_HYPRE_ParaSails(PC pc, PetscViewer viewer)
1130: {
1131: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1132: PetscBool iascii;
1133: const char *symt = 0;
1135: PetscFunctionBegin;
1136: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
1137: if (iascii) {
1138: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE ParaSails preconditioning\n"));
1139: PetscCall(PetscViewerASCIIPrintf(viewer, " nlevels %" PetscInt_FMT "\n", jac->nlevels));
1140: PetscCall(PetscViewerASCIIPrintf(viewer, " threshold %g\n", (double)jac->threshold));
1141: PetscCall(PetscViewerASCIIPrintf(viewer, " filter %g\n", (double)jac->filter));
1142: PetscCall(PetscViewerASCIIPrintf(viewer, " load balance %g\n", (double)jac->loadbal));
1143: PetscCall(PetscViewerASCIIPrintf(viewer, " reuse nonzero structure %s\n", PetscBools[jac->ruse]));
1144: PetscCall(PetscViewerASCIIPrintf(viewer, " print info to screen %s\n", PetscBools[jac->logging]));
1145: if (!jac->symt) symt = "nonsymmetric matrix and preconditioner";
1146: else if (jac->symt == 1) symt = "SPD matrix and preconditioner";
1147: else if (jac->symt == 2) symt = "nonsymmetric matrix but SPD preconditioner";
1148: else SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONG, "Unknown HYPRE ParaSails symmetric option %" PetscInt_FMT, jac->symt);
1149: PetscCall(PetscViewerASCIIPrintf(viewer, " %s\n", symt));
1150: }
1151: PetscFunctionReturn(PETSC_SUCCESS);
1152: }
1154: static PetscErrorCode PCSetFromOptions_HYPRE_AMS(PC pc, PetscOptionItems *PetscOptionsObject)
1155: {
1156: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1157: PetscInt n;
1158: PetscBool flag, flag2, flag3, flag4;
1160: PetscFunctionBegin;
1161: PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE AMS Options");
1162: PetscCall(PetscOptionsInt("-pc_hypre_ams_print_level", "Debugging output level for AMS", "None", jac->as_print, &jac->as_print, &flag));
1163: if (flag) PetscCallExternal(HYPRE_AMSSetPrintLevel, jac->hsolver, jac->as_print);
1164: PetscCall(PetscOptionsInt("-pc_hypre_ams_max_iter", "Maximum number of AMS multigrid iterations within PCApply", "None", jac->as_max_iter, &jac->as_max_iter, &flag));
1165: if (flag) PetscCallExternal(HYPRE_AMSSetMaxIter, jac->hsolver, jac->as_max_iter);
1166: PetscCall(PetscOptionsInt("-pc_hypre_ams_cycle_type", "Cycle type for AMS multigrid", "None", jac->ams_cycle_type, &jac->ams_cycle_type, &flag));
1167: if (flag) PetscCallExternal(HYPRE_AMSSetCycleType, jac->hsolver, jac->ams_cycle_type);
1168: PetscCall(PetscOptionsReal("-pc_hypre_ams_tol", "Error tolerance for AMS multigrid", "None", jac->as_tol, &jac->as_tol, &flag));
1169: if (flag) PetscCallExternal(HYPRE_AMSSetTol, jac->hsolver, jac->as_tol);
1170: PetscCall(PetscOptionsInt("-pc_hypre_ams_relax_type", "Relaxation type for AMS smoother", "None", jac->as_relax_type, &jac->as_relax_type, &flag));
1171: PetscCall(PetscOptionsInt("-pc_hypre_ams_relax_times", "Number of relaxation steps for AMS smoother", "None", jac->as_relax_times, &jac->as_relax_times, &flag2));
1172: PetscCall(PetscOptionsReal("-pc_hypre_ams_relax_weight", "Relaxation weight for AMS smoother", "None", jac->as_relax_weight, &jac->as_relax_weight, &flag3));
1173: PetscCall(PetscOptionsReal("-pc_hypre_ams_omega", "SSOR coefficient for AMS smoother", "None", jac->as_omega, &jac->as_omega, &flag4));
1174: if (flag || flag2 || flag3 || flag4) PetscCallExternal(HYPRE_AMSSetSmoothingOptions, jac->hsolver, jac->as_relax_type, jac->as_relax_times, jac->as_relax_weight, jac->as_omega);
1175: PetscCall(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));
1176: n = 5;
1177: PetscCall(PetscOptionsIntArray("-pc_hypre_ams_amg_alpha_options", "AMG options for vector Poisson", "None", jac->as_amg_alpha_opts, &n, &flag2));
1178: if (flag || flag2) {
1179: PetscCallExternal(HYPRE_AMSSetAlphaAMGOptions, jac->hsolver, jac->as_amg_alpha_opts[0], /* AMG coarsen type */
1180: jac->as_amg_alpha_opts[1], /* AMG agg_levels */
1181: jac->as_amg_alpha_opts[2], /* AMG relax_type */
1182: jac->as_amg_alpha_theta, jac->as_amg_alpha_opts[3], /* AMG interp_type */
1183: jac->as_amg_alpha_opts[4]); /* AMG Pmax */
1184: }
1185: PetscCall(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));
1186: n = 5;
1187: PetscCall(PetscOptionsIntArray("-pc_hypre_ams_amg_beta_options", "AMG options for scalar Poisson solver", "None", jac->as_amg_beta_opts, &n, &flag2));
1188: if (flag || flag2) {
1189: PetscCallExternal(HYPRE_AMSSetBetaAMGOptions, jac->hsolver, jac->as_amg_beta_opts[0], /* AMG coarsen type */
1190: jac->as_amg_beta_opts[1], /* AMG agg_levels */
1191: jac->as_amg_beta_opts[2], /* AMG relax_type */
1192: jac->as_amg_beta_theta, jac->as_amg_beta_opts[3], /* AMG interp_type */
1193: jac->as_amg_beta_opts[4]); /* AMG Pmax */
1194: }
1195: PetscCall(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));
1196: if (flag) { /* override HYPRE's default only if the options is used */
1197: PetscCallExternal(HYPRE_AMSSetProjectionFrequency, jac->hsolver, jac->ams_proj_freq);
1198: }
1199: PetscOptionsHeadEnd();
1200: PetscFunctionReturn(PETSC_SUCCESS);
1201: }
1203: static PetscErrorCode PCView_HYPRE_AMS(PC pc, PetscViewer viewer)
1204: {
1205: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1206: PetscBool iascii;
1208: PetscFunctionBegin;
1209: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
1210: if (iascii) {
1211: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE AMS preconditioning\n"));
1212: PetscCall(PetscViewerASCIIPrintf(viewer, " subspace iterations per application %" PetscInt_FMT "\n", jac->as_max_iter));
1213: PetscCall(PetscViewerASCIIPrintf(viewer, " subspace cycle type %" PetscInt_FMT "\n", jac->ams_cycle_type));
1214: PetscCall(PetscViewerASCIIPrintf(viewer, " subspace iteration tolerance %g\n", (double)jac->as_tol));
1215: PetscCall(PetscViewerASCIIPrintf(viewer, " smoother type %" PetscInt_FMT "\n", jac->as_relax_type));
1216: PetscCall(PetscViewerASCIIPrintf(viewer, " number of smoothing steps %" PetscInt_FMT "\n", jac->as_relax_times));
1217: PetscCall(PetscViewerASCIIPrintf(viewer, " smoother weight %g\n", (double)jac->as_relax_weight));
1218: PetscCall(PetscViewerASCIIPrintf(viewer, " smoother omega %g\n", (double)jac->as_omega));
1219: if (jac->alpha_Poisson) {
1220: PetscCall(PetscViewerASCIIPrintf(viewer, " vector Poisson solver (passed in by user)\n"));
1221: } else {
1222: PetscCall(PetscViewerASCIIPrintf(viewer, " vector Poisson solver (computed) \n"));
1223: }
1224: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG coarsening type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[0]));
1225: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG levels of aggressive coarsening %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[1]));
1226: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG relaxation type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[2]));
1227: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG interpolation type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[3]));
1228: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG max nonzero elements in interpolation rows %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[4]));
1229: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG strength threshold %g\n", (double)jac->as_amg_alpha_theta));
1230: if (!jac->ams_beta_is_zero) {
1231: if (jac->beta_Poisson) {
1232: PetscCall(PetscViewerASCIIPrintf(viewer, " scalar Poisson solver (passed in by user)\n"));
1233: } else {
1234: PetscCall(PetscViewerASCIIPrintf(viewer, " scalar Poisson solver (computed) \n"));
1235: }
1236: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG coarsening type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[0]));
1237: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG levels of aggressive coarsening %" PetscInt_FMT "\n", jac->as_amg_beta_opts[1]));
1238: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG relaxation type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[2]));
1239: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG interpolation type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[3]));
1240: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG max nonzero elements in interpolation rows %" PetscInt_FMT "\n", jac->as_amg_beta_opts[4]));
1241: PetscCall(PetscViewerASCIIPrintf(viewer, " boomerAMG strength threshold %g\n", (double)jac->as_amg_beta_theta));
1242: if (jac->ams_beta_is_zero_part) PetscCall(PetscViewerASCIIPrintf(viewer, " compatible subspace projection frequency %" PetscInt_FMT " (-1 HYPRE uses default)\n", jac->ams_proj_freq));
1243: } else {
1244: PetscCall(PetscViewerASCIIPrintf(viewer, " scalar Poisson solver not used (zero-conductivity everywhere) \n"));
1245: }
1246: }
1247: PetscFunctionReturn(PETSC_SUCCESS);
1248: }
1250: static PetscErrorCode PCSetFromOptions_HYPRE_ADS(PC pc, PetscOptionItems *PetscOptionsObject)
1251: {
1252: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1253: PetscInt n;
1254: PetscBool flag, flag2, flag3, flag4;
1256: PetscFunctionBegin;
1257: PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE ADS Options");
1258: PetscCall(PetscOptionsInt("-pc_hypre_ads_print_level", "Debugging output level for ADS", "None", jac->as_print, &jac->as_print, &flag));
1259: if (flag) PetscCallExternal(HYPRE_ADSSetPrintLevel, jac->hsolver, jac->as_print);
1260: PetscCall(PetscOptionsInt("-pc_hypre_ads_max_iter", "Maximum number of ADS multigrid iterations within PCApply", "None", jac->as_max_iter, &jac->as_max_iter, &flag));
1261: if (flag) PetscCallExternal(HYPRE_ADSSetMaxIter, jac->hsolver, jac->as_max_iter);
1262: PetscCall(PetscOptionsInt("-pc_hypre_ads_cycle_type", "Cycle type for ADS multigrid", "None", jac->ads_cycle_type, &jac->ads_cycle_type, &flag));
1263: if (flag) PetscCallExternal(HYPRE_ADSSetCycleType, jac->hsolver, jac->ads_cycle_type);
1264: PetscCall(PetscOptionsReal("-pc_hypre_ads_tol", "Error tolerance for ADS multigrid", "None", jac->as_tol, &jac->as_tol, &flag));
1265: if (flag) PetscCallExternal(HYPRE_ADSSetTol, jac->hsolver, jac->as_tol);
1266: PetscCall(PetscOptionsInt("-pc_hypre_ads_relax_type", "Relaxation type for ADS smoother", "None", jac->as_relax_type, &jac->as_relax_type, &flag));
1267: PetscCall(PetscOptionsInt("-pc_hypre_ads_relax_times", "Number of relaxation steps for ADS smoother", "None", jac->as_relax_times, &jac->as_relax_times, &flag2));
1268: PetscCall(PetscOptionsReal("-pc_hypre_ads_relax_weight", "Relaxation weight for ADS smoother", "None", jac->as_relax_weight, &jac->as_relax_weight, &flag3));
1269: PetscCall(PetscOptionsReal("-pc_hypre_ads_omega", "SSOR coefficient for ADS smoother", "None", jac->as_omega, &jac->as_omega, &flag4));
1270: if (flag || flag2 || flag3 || flag4) PetscCallExternal(HYPRE_ADSSetSmoothingOptions, jac->hsolver, jac->as_relax_type, jac->as_relax_times, jac->as_relax_weight, jac->as_omega);
1271: PetscCall(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));
1272: n = 5;
1273: PetscCall(PetscOptionsIntArray("-pc_hypre_ads_ams_options", "AMG options for AMS solver inside ADS", "None", jac->as_amg_alpha_opts, &n, &flag2));
1274: PetscCall(PetscOptionsInt("-pc_hypre_ads_ams_cycle_type", "Cycle type for AMS solver inside ADS", "None", jac->ams_cycle_type, &jac->ams_cycle_type, &flag3));
1275: if (flag || flag2 || flag3) {
1276: PetscCallExternal(HYPRE_ADSSetAMSOptions, jac->hsolver, jac->ams_cycle_type, /* AMS cycle type */
1277: jac->as_amg_alpha_opts[0], /* AMG coarsen type */
1278: jac->as_amg_alpha_opts[1], /* AMG agg_levels */
1279: jac->as_amg_alpha_opts[2], /* AMG relax_type */
1280: jac->as_amg_alpha_theta, jac->as_amg_alpha_opts[3], /* AMG interp_type */
1281: jac->as_amg_alpha_opts[4]); /* AMG Pmax */
1282: }
1283: PetscCall(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));
1284: n = 5;
1285: PetscCall(PetscOptionsIntArray("-pc_hypre_ads_amg_options", "AMG options for vector AMG solver inside ADS", "None", jac->as_amg_beta_opts, &n, &flag2));
1286: if (flag || flag2) {
1287: PetscCallExternal(HYPRE_ADSSetAMGOptions, jac->hsolver, jac->as_amg_beta_opts[0], /* AMG coarsen type */
1288: jac->as_amg_beta_opts[1], /* AMG agg_levels */
1289: jac->as_amg_beta_opts[2], /* AMG relax_type */
1290: jac->as_amg_beta_theta, jac->as_amg_beta_opts[3], /* AMG interp_type */
1291: jac->as_amg_beta_opts[4]); /* AMG Pmax */
1292: }
1293: PetscOptionsHeadEnd();
1294: PetscFunctionReturn(PETSC_SUCCESS);
1295: }
1297: static PetscErrorCode PCView_HYPRE_ADS(PC pc, PetscViewer viewer)
1298: {
1299: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1300: PetscBool iascii;
1302: PetscFunctionBegin;
1303: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
1304: if (iascii) {
1305: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE ADS preconditioning\n"));
1306: PetscCall(PetscViewerASCIIPrintf(viewer, " subspace iterations per application %" PetscInt_FMT "\n", jac->as_max_iter));
1307: PetscCall(PetscViewerASCIIPrintf(viewer, " subspace cycle type %" PetscInt_FMT "\n", jac->ads_cycle_type));
1308: PetscCall(PetscViewerASCIIPrintf(viewer, " subspace iteration tolerance %g\n", (double)jac->as_tol));
1309: PetscCall(PetscViewerASCIIPrintf(viewer, " smoother type %" PetscInt_FMT "\n", jac->as_relax_type));
1310: PetscCall(PetscViewerASCIIPrintf(viewer, " number of smoothing steps %" PetscInt_FMT "\n", jac->as_relax_times));
1311: PetscCall(PetscViewerASCIIPrintf(viewer, " smoother weight %g\n", (double)jac->as_relax_weight));
1312: PetscCall(PetscViewerASCIIPrintf(viewer, " smoother omega %g\n", (double)jac->as_omega));
1313: PetscCall(PetscViewerASCIIPrintf(viewer, " AMS solver using boomerAMG\n"));
1314: PetscCall(PetscViewerASCIIPrintf(viewer, " subspace cycle type %" PetscInt_FMT "\n", jac->ams_cycle_type));
1315: PetscCall(PetscViewerASCIIPrintf(viewer, " coarsening type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[0]));
1316: PetscCall(PetscViewerASCIIPrintf(viewer, " levels of aggressive coarsening %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[1]));
1317: PetscCall(PetscViewerASCIIPrintf(viewer, " relaxation type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[2]));
1318: PetscCall(PetscViewerASCIIPrintf(viewer, " interpolation type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[3]));
1319: PetscCall(PetscViewerASCIIPrintf(viewer, " max nonzero elements in interpolation rows %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[4]));
1320: PetscCall(PetscViewerASCIIPrintf(viewer, " strength threshold %g\n", (double)jac->as_amg_alpha_theta));
1321: PetscCall(PetscViewerASCIIPrintf(viewer, " vector Poisson solver using boomerAMG\n"));
1322: PetscCall(PetscViewerASCIIPrintf(viewer, " coarsening type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[0]));
1323: PetscCall(PetscViewerASCIIPrintf(viewer, " levels of aggressive coarsening %" PetscInt_FMT "\n", jac->as_amg_beta_opts[1]));
1324: PetscCall(PetscViewerASCIIPrintf(viewer, " relaxation type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[2]));
1325: PetscCall(PetscViewerASCIIPrintf(viewer, " interpolation type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[3]));
1326: PetscCall(PetscViewerASCIIPrintf(viewer, " max nonzero elements in interpolation rows %" PetscInt_FMT "\n", jac->as_amg_beta_opts[4]));
1327: PetscCall(PetscViewerASCIIPrintf(viewer, " strength threshold %g\n", (double)jac->as_amg_beta_theta));
1328: }
1329: PetscFunctionReturn(PETSC_SUCCESS);
1330: }
1332: static PetscErrorCode PCHYPRESetDiscreteGradient_HYPRE(PC pc, Mat G)
1333: {
1334: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1335: PetscBool ishypre;
1337: PetscFunctionBegin;
1338: PetscCall(PetscObjectTypeCompare((PetscObject)G, MATHYPRE, &ishypre));
1339: if (ishypre) {
1340: PetscCall(PetscObjectReference((PetscObject)G));
1341: PetscCall(MatDestroy(&jac->G));
1342: jac->G = G;
1343: } else {
1344: PetscCall(MatDestroy(&jac->G));
1345: PetscCall(MatConvert(G, MATHYPRE, MAT_INITIAL_MATRIX, &jac->G));
1346: }
1347: PetscFunctionReturn(PETSC_SUCCESS);
1348: }
1350: /*@
1351: PCHYPRESetDiscreteGradient - Set discrete gradient matrix for `PCHYPRE` type of ams or ads
1353: Collective
1355: Input Parameters:
1356: + pc - the preconditioning context
1357: - G - the discrete gradient
1359: Level: intermediate
1361: Notes:
1362: G should have as many rows as the number of edges and as many columns as the number of vertices in the mesh
1364: 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
1366: Developer Note:
1367: This automatically converts the matrix to `MATHYPRE` if it is not already of that type
1369: .seealso: `PCHYPRE`, `PCHYPRESetDiscreteCurl()`
1370: @*/
1371: PetscErrorCode PCHYPRESetDiscreteGradient(PC pc, Mat G)
1372: {
1373: PetscFunctionBegin;
1376: PetscCheckSameComm(pc, 1, G, 2);
1377: PetscTryMethod(pc, "PCHYPRESetDiscreteGradient_C", (PC, Mat), (pc, G));
1378: PetscFunctionReturn(PETSC_SUCCESS);
1379: }
1381: static PetscErrorCode PCHYPRESetDiscreteCurl_HYPRE(PC pc, Mat C)
1382: {
1383: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1384: PetscBool ishypre;
1386: PetscFunctionBegin;
1387: PetscCall(PetscObjectTypeCompare((PetscObject)C, MATHYPRE, &ishypre));
1388: if (ishypre) {
1389: PetscCall(PetscObjectReference((PetscObject)C));
1390: PetscCall(MatDestroy(&jac->C));
1391: jac->C = C;
1392: } else {
1393: PetscCall(MatDestroy(&jac->C));
1394: PetscCall(MatConvert(C, MATHYPRE, MAT_INITIAL_MATRIX, &jac->C));
1395: }
1396: PetscFunctionReturn(PETSC_SUCCESS);
1397: }
1399: /*@
1400: PCHYPRESetDiscreteCurl - Set discrete curl matrx for `PCHYPRE` type of ads
1402: Collective
1404: Input Parameters:
1405: + pc - the preconditioning context
1406: - C - the discrete curl
1408: Level: intermediate
1410: Notes:
1411: C should have as many rows as the number of faces and as many columns as the number of edges in the mesh
1413: 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
1415: Developer Note:
1416: This automatically converts the matrix to `MATHYPRE` if it is not already of that type
1418: If this is only for `PCHYPRE` type of ads it should be called `PCHYPREADSSetDiscreteCurl()`
1420: .seealso: `PCHYPRE`, `PCHYPRESetDiscreteGradient()`
1421: @*/
1422: PetscErrorCode PCHYPRESetDiscreteCurl(PC pc, Mat C)
1423: {
1424: PetscFunctionBegin;
1427: PetscCheckSameComm(pc, 1, C, 2);
1428: PetscTryMethod(pc, "PCHYPRESetDiscreteCurl_C", (PC, Mat), (pc, C));
1429: PetscFunctionReturn(PETSC_SUCCESS);
1430: }
1432: static PetscErrorCode PCHYPRESetInterpolations_HYPRE(PC pc, PetscInt dim, Mat RT_PiFull, Mat RT_Pi[], Mat ND_PiFull, Mat ND_Pi[])
1433: {
1434: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1435: PetscBool ishypre;
1436: PetscInt i;
1437: PetscFunctionBegin;
1439: PetscCall(MatDestroy(&jac->RT_PiFull));
1440: PetscCall(MatDestroy(&jac->ND_PiFull));
1441: for (i = 0; i < 3; ++i) {
1442: PetscCall(MatDestroy(&jac->RT_Pi[i]));
1443: PetscCall(MatDestroy(&jac->ND_Pi[i]));
1444: }
1446: jac->dim = dim;
1447: if (RT_PiFull) {
1448: PetscCall(PetscObjectTypeCompare((PetscObject)RT_PiFull, MATHYPRE, &ishypre));
1449: if (ishypre) {
1450: PetscCall(PetscObjectReference((PetscObject)RT_PiFull));
1451: jac->RT_PiFull = RT_PiFull;
1452: } else {
1453: PetscCall(MatConvert(RT_PiFull, MATHYPRE, MAT_INITIAL_MATRIX, &jac->RT_PiFull));
1454: }
1455: }
1456: if (RT_Pi) {
1457: for (i = 0; i < dim; ++i) {
1458: if (RT_Pi[i]) {
1459: PetscCall(PetscObjectTypeCompare((PetscObject)RT_Pi[i], MATHYPRE, &ishypre));
1460: if (ishypre) {
1461: PetscCall(PetscObjectReference((PetscObject)RT_Pi[i]));
1462: jac->RT_Pi[i] = RT_Pi[i];
1463: } else {
1464: PetscCall(MatConvert(RT_Pi[i], MATHYPRE, MAT_INITIAL_MATRIX, &jac->RT_Pi[i]));
1465: }
1466: }
1467: }
1468: }
1469: if (ND_PiFull) {
1470: PetscCall(PetscObjectTypeCompare((PetscObject)ND_PiFull, MATHYPRE, &ishypre));
1471: if (ishypre) {
1472: PetscCall(PetscObjectReference((PetscObject)ND_PiFull));
1473: jac->ND_PiFull = ND_PiFull;
1474: } else {
1475: PetscCall(MatConvert(ND_PiFull, MATHYPRE, MAT_INITIAL_MATRIX, &jac->ND_PiFull));
1476: }
1477: }
1478: if (ND_Pi) {
1479: for (i = 0; i < dim; ++i) {
1480: if (ND_Pi[i]) {
1481: PetscCall(PetscObjectTypeCompare((PetscObject)ND_Pi[i], MATHYPRE, &ishypre));
1482: if (ishypre) {
1483: PetscCall(PetscObjectReference((PetscObject)ND_Pi[i]));
1484: jac->ND_Pi[i] = ND_Pi[i];
1485: } else {
1486: PetscCall(MatConvert(ND_Pi[i], MATHYPRE, MAT_INITIAL_MATRIX, &jac->ND_Pi[i]));
1487: }
1488: }
1489: }
1490: }
1492: PetscFunctionReturn(PETSC_SUCCESS);
1493: }
1495: /*@
1496: PCHYPRESetInterpolations - Set interpolation matrices for `PCHYPRE` type of ams or ads
1498: Collective
1500: Input Parameters:
1501: + pc - the preconditioning context
1502: . dim - the dimension of the problem, only used in AMS
1503: . RT_PiFull - Raviart-Thomas interpolation matrix
1504: . RT_Pi - x/y/z component of Raviart-Thomas interpolation matrix
1505: . ND_PiFull - Nedelec interpolation matrix
1506: - ND_Pi - x/y/z component of Nedelec interpolation matrix
1508: Level: intermediate
1510: Notes:
1511: For AMS, only Nedelec interpolation matrices are needed, the Raviart-Thomas interpolation matrices can be set to NULL.
1513: For ADS, both type of interpolation matrices are needed.
1515: Developer Note:
1516: This automatically converts the matrix to `MATHYPRE` if it is not already of that type
1518: .seealso: `PCHYPRE`
1519: @*/
1520: PetscErrorCode PCHYPRESetInterpolations(PC pc, PetscInt dim, Mat RT_PiFull, Mat RT_Pi[], Mat ND_PiFull, Mat ND_Pi[])
1521: {
1522: PetscInt i;
1524: PetscFunctionBegin;
1526: if (RT_PiFull) {
1528: PetscCheckSameComm(pc, 1, RT_PiFull, 3);
1529: }
1530: if (RT_Pi) {
1532: for (i = 0; i < dim; ++i) {
1533: if (RT_Pi[i]) {
1535: PetscCheckSameComm(pc, 1, RT_Pi[i], 4);
1536: }
1537: }
1538: }
1539: if (ND_PiFull) {
1541: PetscCheckSameComm(pc, 1, ND_PiFull, 5);
1542: }
1543: if (ND_Pi) {
1545: for (i = 0; i < dim; ++i) {
1546: if (ND_Pi[i]) {
1548: PetscCheckSameComm(pc, 1, ND_Pi[i], 6);
1549: }
1550: }
1551: }
1552: PetscTryMethod(pc, "PCHYPRESetInterpolations_C", (PC, PetscInt, Mat, Mat[], Mat, Mat[]), (pc, dim, RT_PiFull, RT_Pi, ND_PiFull, ND_Pi));
1553: PetscFunctionReturn(PETSC_SUCCESS);
1554: }
1556: static PetscErrorCode PCHYPRESetPoissonMatrix_HYPRE(PC pc, Mat A, PetscBool isalpha)
1557: {
1558: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1559: PetscBool ishypre;
1561: PetscFunctionBegin;
1562: PetscCall(PetscObjectTypeCompare((PetscObject)A, MATHYPRE, &ishypre));
1563: if (ishypre) {
1564: if (isalpha) {
1565: PetscCall(PetscObjectReference((PetscObject)A));
1566: PetscCall(MatDestroy(&jac->alpha_Poisson));
1567: jac->alpha_Poisson = A;
1568: } else {
1569: if (A) {
1570: PetscCall(PetscObjectReference((PetscObject)A));
1571: } else {
1572: jac->ams_beta_is_zero = PETSC_TRUE;
1573: }
1574: PetscCall(MatDestroy(&jac->beta_Poisson));
1575: jac->beta_Poisson = A;
1576: }
1577: } else {
1578: if (isalpha) {
1579: PetscCall(MatDestroy(&jac->alpha_Poisson));
1580: PetscCall(MatConvert(A, MATHYPRE, MAT_INITIAL_MATRIX, &jac->alpha_Poisson));
1581: } else {
1582: if (A) {
1583: PetscCall(MatDestroy(&jac->beta_Poisson));
1584: PetscCall(MatConvert(A, MATHYPRE, MAT_INITIAL_MATRIX, &jac->beta_Poisson));
1585: } else {
1586: PetscCall(MatDestroy(&jac->beta_Poisson));
1587: jac->ams_beta_is_zero = PETSC_TRUE;
1588: }
1589: }
1590: }
1591: PetscFunctionReturn(PETSC_SUCCESS);
1592: }
1594: /*@
1595: PCHYPRESetAlphaPoissonMatrix - Set vector Poisson matrix for `PCHYPRE` of type ams
1597: Collective
1599: Input Parameters:
1600: + pc - the preconditioning context
1601: - A - the matrix
1603: Level: intermediate
1605: Note:
1606: A should be obtained by discretizing the vector valued Poisson problem with linear finite elements
1608: Developer Note:
1609: This automatically converts the matrix to `MATHYPRE` if it is not already of that type
1611: If this is only for `PCHYPRE` type of ams it should be called `PCHYPREAMSSetAlphaPoissonMatrix()`
1613: .seealso: `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetBetaPoissonMatrix()`
1614: @*/
1615: PetscErrorCode PCHYPRESetAlphaPoissonMatrix(PC pc, Mat A)
1616: {
1617: PetscFunctionBegin;
1620: PetscCheckSameComm(pc, 1, A, 2);
1621: PetscTryMethod(pc, "PCHYPRESetPoissonMatrix_C", (PC, Mat, PetscBool), (pc, A, PETSC_TRUE));
1622: PetscFunctionReturn(PETSC_SUCCESS);
1623: }
1625: /*@
1626: PCHYPRESetBetaPoissonMatrix - Set Poisson matrix for `PCHYPRE` of type ams
1628: Collective
1630: Input Parameters:
1631: + pc - the preconditioning context
1632: - A - the matrix, or NULL to turn it off
1634: Level: intermediate
1636: Note:
1637: A should be obtained by discretizing the Poisson problem with linear finite elements.
1639: Developer Note:
1640: This automatically converts the matrix to `MATHYPRE` if it is not already of that type
1642: If this is only for `PCHYPRE` type of ams it should be called `PCHYPREAMSPCHYPRESetBetaPoissonMatrix()`
1644: .seealso: `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetAlphaPoissonMatrix()`
1645: @*/
1646: PetscErrorCode PCHYPRESetBetaPoissonMatrix(PC pc, Mat A)
1647: {
1648: PetscFunctionBegin;
1650: if (A) {
1652: PetscCheckSameComm(pc, 1, A, 2);
1653: }
1654: PetscTryMethod(pc, "PCHYPRESetPoissonMatrix_C", (PC, Mat, PetscBool), (pc, A, PETSC_FALSE));
1655: PetscFunctionReturn(PETSC_SUCCESS);
1656: }
1658: static PetscErrorCode PCHYPRESetEdgeConstantVectors_HYPRE(PC pc, Vec ozz, Vec zoz, Vec zzo)
1659: {
1660: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1662: PetscFunctionBegin;
1663: /* throw away any vector if already set */
1664: PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[0]));
1665: PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[1]));
1666: PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[2]));
1667: PetscCall(VecHYPRE_IJVectorCreate(ozz->map, &jac->constants[0]));
1668: PetscCall(VecHYPRE_IJVectorCopy(ozz, jac->constants[0]));
1669: PetscCall(VecHYPRE_IJVectorCreate(zoz->map, &jac->constants[1]));
1670: PetscCall(VecHYPRE_IJVectorCopy(zoz, jac->constants[1]));
1671: jac->dim = 2;
1672: if (zzo) {
1673: PetscCall(VecHYPRE_IJVectorCreate(zzo->map, &jac->constants[2]));
1674: PetscCall(VecHYPRE_IJVectorCopy(zzo, jac->constants[2]));
1675: jac->dim++;
1676: }
1677: PetscFunctionReturn(PETSC_SUCCESS);
1678: }
1680: /*@
1681: PCHYPRESetEdgeConstantVectors - Set the representation of the constant vector fields in the edge element basis for `PCHYPRE` of type ams
1683: Collective
1685: Input Parameters:
1686: + pc - the preconditioning context
1687: . ozz - vector representing (1,0,0) (or (1,0) in 2D)
1688: . zoz - vector representing (0,1,0) (or (0,1) in 2D)
1689: - zzo - vector representing (0,0,1) (use NULL in 2D)
1691: Level: intermediate
1693: Developer Note:
1694: If this is only for `PCHYPRE` type of ams it should be called `PCHYPREAMSSetEdgeConstantVectors()`
1696: .seealso: `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetAlphaPoissonMatrix()`
1697: @*/
1698: PetscErrorCode PCHYPRESetEdgeConstantVectors(PC pc, Vec ozz, Vec zoz, Vec zzo)
1699: {
1700: PetscFunctionBegin;
1705: PetscCheckSameComm(pc, 1, ozz, 2);
1706: PetscCheckSameComm(pc, 1, zoz, 3);
1707: if (zzo) PetscCheckSameComm(pc, 1, zzo, 4);
1708: PetscTryMethod(pc, "PCHYPRESetEdgeConstantVectors_C", (PC, Vec, Vec, Vec), (pc, ozz, zoz, zzo));
1709: PetscFunctionReturn(PETSC_SUCCESS);
1710: }
1712: static PetscErrorCode PCHYPREAMSSetInteriorNodes_HYPRE(PC pc, Vec interior)
1713: {
1714: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1716: PetscFunctionBegin;
1717: PetscCall(VecHYPRE_IJVectorDestroy(&jac->interior));
1718: PetscCall(VecHYPRE_IJVectorCreate(interior->map, &jac->interior));
1719: PetscCall(VecHYPRE_IJVectorCopy(interior, jac->interior));
1720: jac->ams_beta_is_zero_part = PETSC_TRUE;
1721: PetscFunctionReturn(PETSC_SUCCESS);
1722: }
1724: /*@
1725: PCHYPREAMSSetInteriorNodes - Set the list of interior nodes to a zero-conductivity region for `PCHYPRE` of type ams
1727: Collective
1729: Input Parameters:
1730: + pc - the preconditioning context
1731: - interior - vector. node is interior if its entry in the array is 1.0.
1733: Level: intermediate
1735: Note:
1736: This calls `HYPRE_AMSSetInteriorNodes()`
1738: Developer Note:
1739: If this is only for `PCHYPRE` type of ams it should be called `PCHYPREAMSSetInteriorNodes()`
1741: .seealso: `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetAlphaPoissonMatrix()`
1742: @*/
1743: PetscErrorCode PCHYPREAMSSetInteriorNodes(PC pc, Vec interior)
1744: {
1745: PetscFunctionBegin;
1748: PetscCheckSameComm(pc, 1, interior, 2);
1749: PetscTryMethod(pc, "PCHYPREAMSSetInteriorNodes_C", (PC, Vec), (pc, interior));
1750: PetscFunctionReturn(PETSC_SUCCESS);
1751: }
1753: static PetscErrorCode PCSetCoordinates_HYPRE(PC pc, PetscInt dim, PetscInt nloc, PetscReal *coords)
1754: {
1755: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1756: Vec tv;
1757: PetscInt i;
1759: PetscFunctionBegin;
1760: /* throw away any coordinate vector if already set */
1761: PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[0]));
1762: PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[1]));
1763: PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[2]));
1764: jac->dim = dim;
1766: /* compute IJ vector for coordinates */
1767: PetscCall(VecCreate(PetscObjectComm((PetscObject)pc), &tv));
1768: PetscCall(VecSetType(tv, VECSTANDARD));
1769: PetscCall(VecSetSizes(tv, nloc, PETSC_DECIDE));
1770: for (i = 0; i < dim; i++) {
1771: PetscScalar *array;
1772: PetscInt j;
1774: PetscCall(VecHYPRE_IJVectorCreate(tv->map, &jac->coords[i]));
1775: PetscCall(VecGetArrayWrite(tv, &array));
1776: for (j = 0; j < nloc; j++) array[j] = coords[j * dim + i];
1777: PetscCall(VecRestoreArrayWrite(tv, &array));
1778: PetscCall(VecHYPRE_IJVectorCopy(tv, jac->coords[i]));
1779: }
1780: PetscCall(VecDestroy(&tv));
1781: PetscFunctionReturn(PETSC_SUCCESS);
1782: }
1784: static PetscErrorCode PCHYPREGetType_HYPRE(PC pc, const char *name[])
1785: {
1786: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1788: PetscFunctionBegin;
1789: *name = jac->hypre_type;
1790: PetscFunctionReturn(PETSC_SUCCESS);
1791: }
1793: static PetscErrorCode PCHYPRESetType_HYPRE(PC pc, const char name[])
1794: {
1795: PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1796: PetscBool flag;
1798: PetscFunctionBegin;
1799: if (jac->hypre_type) {
1800: PetscCall(PetscStrcmp(jac->hypre_type, name, &flag));
1801: PetscCheck(flag, PetscObjectComm((PetscObject)pc), PETSC_ERR_ORDER, "Cannot reset the HYPRE preconditioner type once it has been set");
1802: PetscFunctionReturn(PETSC_SUCCESS);
1803: } else {
1804: PetscCall(PetscStrallocpy(name, &jac->hypre_type));
1805: }
1807: jac->maxiter = PETSC_DEFAULT;
1808: jac->tol = PETSC_DEFAULT;
1809: jac->printstatistics = PetscLogPrintInfo;
1811: PetscCall(PetscStrcmp("pilut", jac->hypre_type, &flag));
1812: if (flag) {
1813: PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &jac->comm_hypre));
1814: PetscCallExternal(HYPRE_ParCSRPilutCreate, jac->comm_hypre, &jac->hsolver);
1815: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Pilut;
1816: pc->ops->view = PCView_HYPRE_Pilut;
1817: jac->destroy = HYPRE_ParCSRPilutDestroy;
1818: jac->setup = HYPRE_ParCSRPilutSetup;
1819: jac->solve = HYPRE_ParCSRPilutSolve;
1820: jac->factorrowsize = PETSC_DEFAULT;
1821: PetscFunctionReturn(PETSC_SUCCESS);
1822: }
1823: PetscCall(PetscStrcmp("euclid", jac->hypre_type, &flag));
1824: if (flag) {
1825: #if defined(PETSC_USE_64BIT_INDICES)
1826: SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_SUP, "Hypre Euclid does not support 64-bit indices");
1827: #endif
1828: PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &jac->comm_hypre));
1829: PetscCallExternal(HYPRE_EuclidCreate, jac->comm_hypre, &jac->hsolver);
1830: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Euclid;
1831: pc->ops->view = PCView_HYPRE_Euclid;
1832: jac->destroy = HYPRE_EuclidDestroy;
1833: jac->setup = HYPRE_EuclidSetup;
1834: jac->solve = HYPRE_EuclidSolve;
1835: jac->factorrowsize = PETSC_DEFAULT;
1836: jac->eu_level = PETSC_DEFAULT; /* default */
1837: PetscFunctionReturn(PETSC_SUCCESS);
1838: }
1839: PetscCall(PetscStrcmp("parasails", jac->hypre_type, &flag));
1840: if (flag) {
1841: PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &jac->comm_hypre));
1842: PetscCallExternal(HYPRE_ParaSailsCreate, jac->comm_hypre, &jac->hsolver);
1843: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ParaSails;
1844: pc->ops->view = PCView_HYPRE_ParaSails;
1845: jac->destroy = HYPRE_ParaSailsDestroy;
1846: jac->setup = HYPRE_ParaSailsSetup;
1847: jac->solve = HYPRE_ParaSailsSolve;
1848: /* initialize */
1849: jac->nlevels = 1;
1850: jac->threshold = .1;
1851: jac->filter = .1;
1852: jac->loadbal = 0;
1853: if (PetscLogPrintInfo) jac->logging = (int)PETSC_TRUE;
1854: else jac->logging = (int)PETSC_FALSE;
1856: jac->ruse = (int)PETSC_FALSE;
1857: jac->symt = 0;
1858: PetscCallExternal(HYPRE_ParaSailsSetParams, jac->hsolver, jac->threshold, jac->nlevels);
1859: PetscCallExternal(HYPRE_ParaSailsSetFilter, jac->hsolver, jac->filter);
1860: PetscCallExternal(HYPRE_ParaSailsSetLoadbal, jac->hsolver, jac->loadbal);
1861: PetscCallExternal(HYPRE_ParaSailsSetLogging, jac->hsolver, jac->logging);
1862: PetscCallExternal(HYPRE_ParaSailsSetReuse, jac->hsolver, jac->ruse);
1863: PetscCallExternal(HYPRE_ParaSailsSetSym, jac->hsolver, jac->symt);
1864: PetscFunctionReturn(PETSC_SUCCESS);
1865: }
1866: PetscCall(PetscStrcmp("boomeramg", jac->hypre_type, &flag));
1867: if (flag) {
1868: PetscCallExternal(HYPRE_BoomerAMGCreate, &jac->hsolver);
1869: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_BoomerAMG;
1870: pc->ops->view = PCView_HYPRE_BoomerAMG;
1871: pc->ops->applytranspose = PCApplyTranspose_HYPRE_BoomerAMG;
1872: pc->ops->applyrichardson = PCApplyRichardson_HYPRE_BoomerAMG;
1873: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCGetInterpolations_C", PCGetInterpolations_BoomerAMG));
1874: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCGetCoarseOperators_C", PCGetCoarseOperators_BoomerAMG));
1875: jac->destroy = HYPRE_BoomerAMGDestroy;
1876: jac->setup = HYPRE_BoomerAMGSetup;
1877: jac->solve = HYPRE_BoomerAMGSolve;
1878: jac->applyrichardson = PETSC_FALSE;
1879: /* these defaults match the hypre defaults */
1880: jac->cycletype = 1;
1881: jac->maxlevels = 25;
1882: jac->maxiter = 1;
1883: jac->tol = 0.0; /* tolerance of zero indicates use as preconditioner (suppresses convergence errors) */
1884: jac->truncfactor = 0.0;
1885: jac->strongthreshold = .25;
1886: jac->maxrowsum = .9;
1887: jac->coarsentype = 6;
1888: jac->measuretype = 0;
1889: jac->gridsweeps[0] = jac->gridsweeps[1] = jac->gridsweeps[2] = 1;
1890: jac->smoothtype = -1; /* Not set by default */
1891: jac->smoothnumlevels = 25;
1892: jac->eu_level = 0;
1893: jac->eu_droptolerance = 0;
1894: jac->eu_bj = 0;
1895: jac->relaxtype[0] = jac->relaxtype[1] = 6; /* Defaults to SYMMETRIC since in PETSc we are using a PC - most likely with CG */
1896: jac->relaxtype[2] = 9; /*G.E. */
1897: jac->relaxweight = 1.0;
1898: jac->outerrelaxweight = 1.0;
1899: jac->relaxorder = 1;
1900: jac->interptype = 0;
1901: jac->Rtype = 0;
1902: jac->Rstrongthreshold = 0.25;
1903: jac->Rfilterthreshold = 0.0;
1904: jac->Adroptype = -1;
1905: jac->Adroptol = 0.0;
1906: jac->agg_nl = 0;
1907: jac->agg_interptype = 4;
1908: jac->pmax = 0;
1909: jac->truncfactor = 0.0;
1910: jac->agg_num_paths = 1;
1911: jac->maxc = 9;
1912: jac->minc = 1;
1913: jac->nodal_coarsening = 0;
1914: jac->nodal_coarsening_diag = 0;
1915: jac->vec_interp_variant = 0;
1916: jac->vec_interp_qmax = 0;
1917: jac->vec_interp_smooth = PETSC_FALSE;
1918: jac->interp_refine = 0;
1919: jac->nodal_relax = PETSC_FALSE;
1920: jac->nodal_relax_levels = 1;
1921: jac->rap2 = 0;
1923: /* GPU defaults
1924: from https://hypre.readthedocs.io/en/latest/solvers-boomeramg.html#gpu-supported-options
1925: and /src/parcsr_ls/par_amg.c */
1926: #if defined(PETSC_HAVE_HYPRE_DEVICE)
1927: jac->keeptranspose = PETSC_TRUE;
1928: jac->mod_rap2 = 1;
1929: jac->coarsentype = 8;
1930: jac->relaxorder = 0;
1931: jac->interptype = 6;
1932: jac->relaxtype[0] = 18;
1933: jac->relaxtype[1] = 18;
1934: jac->agg_interptype = 7;
1935: #else
1936: jac->keeptranspose = PETSC_FALSE;
1937: jac->mod_rap2 = 0;
1938: #endif
1939: PetscCallExternal(HYPRE_BoomerAMGSetCycleType, jac->hsolver, jac->cycletype);
1940: PetscCallExternal(HYPRE_BoomerAMGSetMaxLevels, jac->hsolver, jac->maxlevels);
1941: PetscCallExternal(HYPRE_BoomerAMGSetMaxIter, jac->hsolver, jac->maxiter);
1942: PetscCallExternal(HYPRE_BoomerAMGSetTol, jac->hsolver, jac->tol);
1943: PetscCallExternal(HYPRE_BoomerAMGSetTruncFactor, jac->hsolver, jac->truncfactor);
1944: PetscCallExternal(HYPRE_BoomerAMGSetStrongThreshold, jac->hsolver, jac->strongthreshold);
1945: PetscCallExternal(HYPRE_BoomerAMGSetMaxRowSum, jac->hsolver, jac->maxrowsum);
1946: PetscCallExternal(HYPRE_BoomerAMGSetCoarsenType, jac->hsolver, jac->coarsentype);
1947: PetscCallExternal(HYPRE_BoomerAMGSetMeasureType, jac->hsolver, jac->measuretype);
1948: PetscCallExternal(HYPRE_BoomerAMGSetRelaxOrder, jac->hsolver, jac->relaxorder);
1949: PetscCallExternal(HYPRE_BoomerAMGSetInterpType, jac->hsolver, jac->interptype);
1950: PetscCallExternal(HYPRE_BoomerAMGSetAggNumLevels, jac->hsolver, jac->agg_nl);
1951: PetscCallExternal(HYPRE_BoomerAMGSetAggInterpType, jac->hsolver, jac->agg_interptype);
1952: PetscCallExternal(HYPRE_BoomerAMGSetPMaxElmts, jac->hsolver, jac->pmax);
1953: PetscCallExternal(HYPRE_BoomerAMGSetNumPaths, jac->hsolver, jac->agg_num_paths);
1954: PetscCallExternal(HYPRE_BoomerAMGSetRelaxType, jac->hsolver, jac->relaxtype[0]); /* defaults coarse to 9 */
1955: PetscCallExternal(HYPRE_BoomerAMGSetNumSweeps, jac->hsolver, jac->gridsweeps[0]); /* defaults coarse to 1 */
1956: PetscCallExternal(HYPRE_BoomerAMGSetMaxCoarseSize, jac->hsolver, jac->maxc);
1957: PetscCallExternal(HYPRE_BoomerAMGSetMinCoarseSize, jac->hsolver, jac->minc);
1958: /* GPU */
1959: #if PETSC_PKG_HYPRE_VERSION_GE(2, 18, 0)
1960: PetscCallExternal(HYPRE_BoomerAMGSetKeepTranspose, jac->hsolver, jac->keeptranspose ? 1 : 0);
1961: PetscCallExternal(HYPRE_BoomerAMGSetRAP2, jac->hsolver, jac->rap2);
1962: PetscCallExternal(HYPRE_BoomerAMGSetModuleRAP2, jac->hsolver, jac->mod_rap2);
1963: #endif
1965: /* AIR */
1966: #if PETSC_PKG_HYPRE_VERSION_GE(2, 18, 0)
1967: PetscCallExternal(HYPRE_BoomerAMGSetRestriction, jac->hsolver, jac->Rtype);
1968: PetscCallExternal(HYPRE_BoomerAMGSetStrongThresholdR, jac->hsolver, jac->Rstrongthreshold);
1969: PetscCallExternal(HYPRE_BoomerAMGSetFilterThresholdR, jac->hsolver, jac->Rfilterthreshold);
1970: PetscCallExternal(HYPRE_BoomerAMGSetADropTol, jac->hsolver, jac->Adroptol);
1971: PetscCallExternal(HYPRE_BoomerAMGSetADropType, jac->hsolver, jac->Adroptype);
1972: #endif
1973: PetscFunctionReturn(PETSC_SUCCESS);
1974: }
1975: PetscCall(PetscStrcmp("ams", jac->hypre_type, &flag));
1976: if (flag) {
1977: PetscCallExternal(HYPRE_AMSCreate, &jac->hsolver);
1978: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_AMS;
1979: pc->ops->view = PCView_HYPRE_AMS;
1980: jac->destroy = HYPRE_AMSDestroy;
1981: jac->setup = HYPRE_AMSSetup;
1982: jac->solve = HYPRE_AMSSolve;
1983: jac->coords[0] = NULL;
1984: jac->coords[1] = NULL;
1985: jac->coords[2] = NULL;
1986: jac->interior = NULL;
1987: /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE AMS */
1988: jac->as_print = 0;
1989: jac->as_max_iter = 1; /* used as a preconditioner */
1990: jac->as_tol = 0.; /* used as a preconditioner */
1991: jac->ams_cycle_type = 13;
1992: /* Smoothing options */
1993: jac->as_relax_type = 2;
1994: jac->as_relax_times = 1;
1995: jac->as_relax_weight = 1.0;
1996: jac->as_omega = 1.0;
1997: /* Vector valued Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
1998: jac->as_amg_alpha_opts[0] = 10;
1999: jac->as_amg_alpha_opts[1] = 1;
2000: jac->as_amg_alpha_opts[2] = 6;
2001: jac->as_amg_alpha_opts[3] = 6;
2002: jac->as_amg_alpha_opts[4] = 4;
2003: jac->as_amg_alpha_theta = 0.25;
2004: /* Scalar Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
2005: jac->as_amg_beta_opts[0] = 10;
2006: jac->as_amg_beta_opts[1] = 1;
2007: jac->as_amg_beta_opts[2] = 6;
2008: jac->as_amg_beta_opts[3] = 6;
2009: jac->as_amg_beta_opts[4] = 4;
2010: jac->as_amg_beta_theta = 0.25;
2011: PetscCallExternal(HYPRE_AMSSetPrintLevel, jac->hsolver, jac->as_print);
2012: PetscCallExternal(HYPRE_AMSSetMaxIter, jac->hsolver, jac->as_max_iter);
2013: PetscCallExternal(HYPRE_AMSSetCycleType, jac->hsolver, jac->ams_cycle_type);
2014: PetscCallExternal(HYPRE_AMSSetTol, jac->hsolver, jac->as_tol);
2015: PetscCallExternal(HYPRE_AMSSetSmoothingOptions, jac->hsolver, jac->as_relax_type, jac->as_relax_times, jac->as_relax_weight, jac->as_omega);
2016: PetscCallExternal(HYPRE_AMSSetAlphaAMGOptions, jac->hsolver, jac->as_amg_alpha_opts[0], /* AMG coarsen type */
2017: jac->as_amg_alpha_opts[1], /* AMG agg_levels */
2018: jac->as_amg_alpha_opts[2], /* AMG relax_type */
2019: jac->as_amg_alpha_theta, jac->as_amg_alpha_opts[3], /* AMG interp_type */
2020: jac->as_amg_alpha_opts[4]); /* AMG Pmax */
2021: PetscCallExternal(HYPRE_AMSSetBetaAMGOptions, jac->hsolver, jac->as_amg_beta_opts[0], /* AMG coarsen type */
2022: jac->as_amg_beta_opts[1], /* AMG agg_levels */
2023: jac->as_amg_beta_opts[2], /* AMG relax_type */
2024: jac->as_amg_beta_theta, jac->as_amg_beta_opts[3], /* AMG interp_type */
2025: jac->as_amg_beta_opts[4]); /* AMG Pmax */
2026: /* Zero conductivity */
2027: jac->ams_beta_is_zero = PETSC_FALSE;
2028: jac->ams_beta_is_zero_part = PETSC_FALSE;
2029: PetscFunctionReturn(PETSC_SUCCESS);
2030: }
2031: PetscCall(PetscStrcmp("ads", jac->hypre_type, &flag));
2032: if (flag) {
2033: PetscCallExternal(HYPRE_ADSCreate, &jac->hsolver);
2034: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ADS;
2035: pc->ops->view = PCView_HYPRE_ADS;
2036: jac->destroy = HYPRE_ADSDestroy;
2037: jac->setup = HYPRE_ADSSetup;
2038: jac->solve = HYPRE_ADSSolve;
2039: jac->coords[0] = NULL;
2040: jac->coords[1] = NULL;
2041: jac->coords[2] = NULL;
2042: /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE ADS */
2043: jac->as_print = 0;
2044: jac->as_max_iter = 1; /* used as a preconditioner */
2045: jac->as_tol = 0.; /* used as a preconditioner */
2046: jac->ads_cycle_type = 13;
2047: /* Smoothing options */
2048: jac->as_relax_type = 2;
2049: jac->as_relax_times = 1;
2050: jac->as_relax_weight = 1.0;
2051: jac->as_omega = 1.0;
2052: /* AMS solver parameters: cycle_type, coarsen type, agg_levels, relax_type, interp_type, Pmax */
2053: jac->ams_cycle_type = 14;
2054: jac->as_amg_alpha_opts[0] = 10;
2055: jac->as_amg_alpha_opts[1] = 1;
2056: jac->as_amg_alpha_opts[2] = 6;
2057: jac->as_amg_alpha_opts[3] = 6;
2058: jac->as_amg_alpha_opts[4] = 4;
2059: jac->as_amg_alpha_theta = 0.25;
2060: /* Vector Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
2061: jac->as_amg_beta_opts[0] = 10;
2062: jac->as_amg_beta_opts[1] = 1;
2063: jac->as_amg_beta_opts[2] = 6;
2064: jac->as_amg_beta_opts[3] = 6;
2065: jac->as_amg_beta_opts[4] = 4;
2066: jac->as_amg_beta_theta = 0.25;
2067: PetscCallExternal(HYPRE_ADSSetPrintLevel, jac->hsolver, jac->as_print);
2068: PetscCallExternal(HYPRE_ADSSetMaxIter, jac->hsolver, jac->as_max_iter);
2069: PetscCallExternal(HYPRE_ADSSetCycleType, jac->hsolver, jac->ams_cycle_type);
2070: PetscCallExternal(HYPRE_ADSSetTol, jac->hsolver, jac->as_tol);
2071: PetscCallExternal(HYPRE_ADSSetSmoothingOptions, jac->hsolver, jac->as_relax_type, jac->as_relax_times, jac->as_relax_weight, jac->as_omega);
2072: PetscCallExternal(HYPRE_ADSSetAMSOptions, jac->hsolver, jac->ams_cycle_type, /* AMG coarsen type */
2073: jac->as_amg_alpha_opts[0], /* AMG coarsen type */
2074: jac->as_amg_alpha_opts[1], /* AMG agg_levels */
2075: jac->as_amg_alpha_opts[2], /* AMG relax_type */
2076: jac->as_amg_alpha_theta, jac->as_amg_alpha_opts[3], /* AMG interp_type */
2077: jac->as_amg_alpha_opts[4]); /* AMG Pmax */
2078: PetscCallExternal(HYPRE_ADSSetAMGOptions, jac->hsolver, jac->as_amg_beta_opts[0], /* AMG coarsen type */
2079: jac->as_amg_beta_opts[1], /* AMG agg_levels */
2080: jac->as_amg_beta_opts[2], /* AMG relax_type */
2081: jac->as_amg_beta_theta, jac->as_amg_beta_opts[3], /* AMG interp_type */
2082: jac->as_amg_beta_opts[4]); /* AMG Pmax */
2083: PetscFunctionReturn(PETSC_SUCCESS);
2084: }
2085: PetscCall(PetscFree(jac->hypre_type));
2087: jac->hypre_type = NULL;
2088: SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown HYPRE preconditioner %s; Choices are euclid, pilut, parasails, boomeramg, ams", name);
2089: }
2091: /*
2092: It only gets here if the HYPRE type has not been set before the call to
2093: ...SetFromOptions() which actually is most of the time
2094: */
2095: PetscErrorCode PCSetFromOptions_HYPRE(PC pc, PetscOptionItems *PetscOptionsObject)
2096: {
2097: PetscInt indx;
2098: const char *type[] = {"euclid", "pilut", "parasails", "boomeramg", "ams", "ads"};
2099: PetscBool flg;
2101: PetscFunctionBegin;
2102: PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE preconditioner options");
2103: PetscCall(PetscOptionsEList("-pc_hypre_type", "HYPRE preconditioner type", "PCHYPRESetType", type, PETSC_STATIC_ARRAY_LENGTH(type), "boomeramg", &indx, &flg));
2104: if (flg) {
2105: PetscCall(PCHYPRESetType_HYPRE(pc, type[indx]));
2106: } else {
2107: PetscCall(PCHYPRESetType_HYPRE(pc, "boomeramg"));
2108: }
2109: PetscTryTypeMethod(pc, setfromoptions, PetscOptionsObject);
2110: PetscOptionsHeadEnd();
2111: PetscFunctionReturn(PETSC_SUCCESS);
2112: }
2114: /*@C
2115: PCHYPRESetType - Sets which hypre preconditioner you wish to use
2117: Input Parameters:
2118: + pc - the preconditioner context
2119: - name - either euclid, pilut, parasails, boomeramg, ams, ads
2121: Options Database Key:
2122: -pc_hypre_type - One of euclid, pilut, parasails, boomeramg, ams, ads
2124: Level: intermediate
2126: .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`, `PCHYPRE`
2127: @*/
2128: PetscErrorCode PCHYPRESetType(PC pc, const char name[])
2129: {
2130: PetscFunctionBegin;
2133: PetscTryMethod(pc, "PCHYPRESetType_C", (PC, const char[]), (pc, name));
2134: PetscFunctionReturn(PETSC_SUCCESS);
2135: }
2137: /*@C
2138: PCHYPREGetType - Gets which hypre preconditioner you are using
2140: Input Parameter:
2141: . pc - the preconditioner context
2143: Output Parameter:
2144: . name - either euclid, pilut, parasails, boomeramg, ams, ads
2146: Level: intermediate
2148: .seealso: `PCCreate()`, `PCHYPRESetType()`, `PCType`, `PC`, `PCHYPRE`
2149: @*/
2150: PetscErrorCode PCHYPREGetType(PC pc, const char *name[])
2151: {
2152: PetscFunctionBegin;
2155: PetscTryMethod(pc, "PCHYPREGetType_C", (PC, const char *[]), (pc, name));
2156: PetscFunctionReturn(PETSC_SUCCESS);
2157: }
2159: /*@C
2160: PCMGGalerkinSetMatProductAlgorithm - Set type of SpGEMM for hypre to use on GPUs
2162: Logically Collective
2164: Input Parameters:
2165: + pc - the hypre context
2166: - type - one of 'cusparse', 'hypre'
2168: Options Database Key:
2169: . -pc_mg_galerkin_mat_product_algorithm <cusparse,hypre> - Type of SpGEMM to use in hypre
2171: Level: intermediate
2173: Developer Note:
2174: How the name starts with `PCMG`, should it not be `PCHYPREBoomerAMG`?
2176: .seealso: `PCHYPRE`, `PCMGGalerkinGetMatProductAlgorithm()`
2177: @*/
2178: PetscErrorCode PCMGGalerkinSetMatProductAlgorithm(PC pc, const char name[])
2179: {
2180: PetscFunctionBegin;
2182: PetscTryMethod(pc, "PCMGGalerkinSetMatProductAlgorithm_C", (PC, const char[]), (pc, name));
2183: PetscFunctionReturn(PETSC_SUCCESS);
2184: }
2186: /*@C
2187: PCMGGalerkinGetMatProductAlgorithm - Get type of SpGEMM for hypre to use on GPUs
2189: Not Collective
2191: Input Parameter:
2192: . pc - the multigrid context
2194: Output Parameter:
2195: . name - one of 'cusparse', 'hypre'
2197: Level: intermediate
2199: .seealso: `PCHYPRE`, ``PCMGGalerkinSetMatProductAlgorithm()`
2200: @*/
2201: PetscErrorCode PCMGGalerkinGetMatProductAlgorithm(PC pc, const char *name[])
2202: {
2203: PetscFunctionBegin;
2205: PetscTryMethod(pc, "PCMGGalerkinGetMatProductAlgorithm_C", (PC, const char *[]), (pc, name));
2206: PetscFunctionReturn(PETSC_SUCCESS);
2207: }
2209: /*MC
2210: PCHYPRE - Allows you to use the matrix element based preconditioners in the LLNL package hypre as PETSc `PC`
2212: Options Database Keys:
2213: + -pc_hypre_type - One of `euclid`, `pilut`, `parasails`, `boomeramg`, `ams`, or `ads`
2214: . -pc_hypre_boomeramg_nodal_coarsen <n> - where n is from 1 to 6 (see `HYPRE_BOOMERAMGSetNodal()`)
2215: . -pc_hypre_boomeramg_vec_interp_variant <v> - where v is from 1 to 3 (see `HYPRE_BoomerAMGSetInterpVecVariant()`)
2216: - Many others, run with `-pc_type hypre` `-pc_hypre_type XXX` `-help` to see options for the XXX preconditioner
2218: Level: intermediate
2220: Notes:
2221: Apart from `-pc_hypre_type` (for which there is `PCHYPRESetType()`),
2222: the many hypre options can ONLY be set via the options database (e.g. the command line
2223: or with `PetscOptionsSetValue()`, there are no functions to set them)
2225: The options `-pc_hypre_boomeramg_max_iter` and `-pc_hypre_boomeramg_tol` refer to the number of iterations
2226: (V-cycles) and tolerance that boomerAMG does EACH time it is called. So for example, if
2227: `-pc_hypre_boomeramg_max_iter` is set to 2 then 2-V-cycles are being used to define the preconditioner
2228: (`-pc_hypre_boomeramg_tol` should be set to 0.0 - the default - to strictly use a fixed number of
2229: iterations per hypre call). `-ksp_max_it` and `-ksp_rtol` STILL determine the total number of iterations
2230: and tolerance for the Krylov solver. For example, if `-pc_hypre_boomeramg_max_iter` is 2 and `-ksp_max_it` is 10
2231: then AT MOST twenty V-cycles of boomeramg will be used.
2233: Note that the option `-pc_hypre_boomeramg_relax_type_all` defaults to symmetric relaxation
2234: (symmetric-SOR/Jacobi), which is required for Krylov solvers like CG that expect symmetry.
2235: Otherwise, you may want to use `-pc_hypre_boomeramg_relax_type_all SOR/Jacobi`.
2237: `MatSetNearNullSpace()` - if you provide a near null space to your matrix it is ignored by hypre UNLESS you also use
2238: the following two options: `-pc_hypre_boomeramg_nodal_coarsen <n> -pc_hypre_boomeramg_vec_interp_variant <v>`
2240: See `PCPFMG`, `PCSMG`, and `PCSYSPFMG` for access to hypre's other (nonalgebraic) multigrid solvers
2242: For `PCHYPRE` type of `ams` or `ads` auxiliary data must be provided to the preconditioner with `PCHYPRESetDiscreteGradient()`,
2243: `PCHYPRESetDiscreteCurl()`, `PCHYPRESetInterpolations()`, `PCHYPRESetAlphaPoissonMatrix()`, `PCHYPRESetBetaPoissonMatrix()`, `PCHYPRESetEdgeConstantVectors()`,
2244: `PCHYPREAMSSetInteriorNodes()`
2246: Sometimes people want to try algebraic multigrid as a "standalone" solver, that is not accelerating it with a Krylov method. Though we generally do not recommend this
2247: since it is usually slower, one should use a `KSPType` of `KSPRICHARDSON`
2248: (or equivalently `-ksp_type richardson`) to achieve this. Using `KSPPREONLY` will not work since it only applies a single cycle of multigrid.
2250: PETSc provides its own geometric and algebraic multigrid solvers `PCMG` and `PCGAMG`, also see `PCHMG` which is useful for certain multicomponent problems
2252: GPU Notes:
2253: To configure hypre BoomerAMG so that it can utilize NVIDIA GPUs run ./configure --download-hypre --with-cuda
2254: Then pass `VECCUDA` vectors and `MATAIJCUSPARSE` matrices to the solvers and PETSc will automatically utilize hypre's GPU solvers.
2256: To configure hypre BoomerAMG so that it can utilize AMD GPUs run ./configure --download-hypre --with-hip
2257: Then pass `VECHIP` vectors to the solvers and PETSc will automatically utilize hypre's GPU solvers.
2259: .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`, `PCHYPRESetType()`, `PCPFMG`, `PCGAMG`, `PCSYSPFMG`, `PCSMG`, `PCHYPRESetDiscreteGradient()`,
2260: `PCHYPRESetDiscreteCurl()`, `PCHYPRESetInterpolations()`, `PCHYPRESetAlphaPoissonMatrix()`, `PCHYPRESetBetaPoissonMatrix()`, `PCHYPRESetEdgeConstantVectors()`,
2261: PCHYPREAMSSetInteriorNodes()
2262: M*/
2264: PETSC_EXTERN PetscErrorCode PCCreate_HYPRE(PC pc)
2265: {
2266: PC_HYPRE *jac;
2268: PetscFunctionBegin;
2269: PetscCall(PetscNew(&jac));
2271: pc->data = jac;
2272: pc->ops->reset = PCReset_HYPRE;
2273: pc->ops->destroy = PCDestroy_HYPRE;
2274: pc->ops->setfromoptions = PCSetFromOptions_HYPRE;
2275: pc->ops->setup = PCSetUp_HYPRE;
2276: pc->ops->apply = PCApply_HYPRE;
2277: jac->comm_hypre = MPI_COMM_NULL;
2278: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetType_C", PCHYPRESetType_HYPRE));
2279: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPREGetType_C", PCHYPREGetType_HYPRE));
2280: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSetCoordinates_C", PCSetCoordinates_HYPRE));
2281: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetDiscreteGradient_C", PCHYPRESetDiscreteGradient_HYPRE));
2282: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetDiscreteCurl_C", PCHYPRESetDiscreteCurl_HYPRE));
2283: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetInterpolations_C", PCHYPRESetInterpolations_HYPRE));
2284: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetEdgeConstantVectors_C", PCHYPRESetEdgeConstantVectors_HYPRE));
2285: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPREAMSSetInteriorNodes_C", PCHYPREAMSSetInteriorNodes_HYPRE));
2286: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetPoissonMatrix_C", PCHYPRESetPoissonMatrix_HYPRE));
2287: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCMGGalerkinSetMatProductAlgorithm_C", PCMGGalerkinSetMatProductAlgorithm_HYPRE_BoomerAMG));
2288: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCMGGalerkinGetMatProductAlgorithm_C", PCMGGalerkinGetMatProductAlgorithm_HYPRE_BoomerAMG));
2289: #if defined(PETSC_HAVE_HYPRE_DEVICE)
2290: #if defined(HYPRE_USING_HIP)
2291: PetscCall(PetscDeviceInitialize(PETSC_DEVICE_HIP));
2292: #endif
2293: #if defined(HYPRE_USING_CUDA)
2294: PetscCall(PetscDeviceInitialize(PETSC_DEVICE_CUDA));
2295: #endif
2296: #endif
2297: PetscFunctionReturn(PETSC_SUCCESS);
2298: }
2300: typedef struct {
2301: MPI_Comm hcomm; /* does not share comm with HYPRE_StructMatrix because need to create solver before getting matrix */
2302: HYPRE_StructSolver hsolver;
2304: /* keep copy of PFMG options used so may view them */
2305: PetscInt its;
2306: double tol;
2307: PetscInt relax_type;
2308: PetscInt rap_type;
2309: PetscInt num_pre_relax, num_post_relax;
2310: PetscInt max_levels;
2311: PetscInt skip_relax;
2312: PetscBool print_statistics;
2313: } PC_PFMG;
2315: PetscErrorCode PCDestroy_PFMG(PC pc)
2316: {
2317: PC_PFMG *ex = (PC_PFMG *)pc->data;
2319: PetscFunctionBegin;
2320: if (ex->hsolver) PetscCallExternal(HYPRE_StructPFMGDestroy, ex->hsolver);
2321: PetscCall(PetscCommRestoreComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
2322: PetscCall(PetscFree(pc->data));
2323: PetscFunctionReturn(PETSC_SUCCESS);
2324: }
2326: static const char *PFMGRelaxType[] = {"Jacobi", "Weighted-Jacobi", "symmetric-Red/Black-Gauss-Seidel", "Red/Black-Gauss-Seidel"};
2327: static const char *PFMGRAPType[] = {"Galerkin", "non-Galerkin"};
2329: PetscErrorCode PCView_PFMG(PC pc, PetscViewer viewer)
2330: {
2331: PetscBool iascii;
2332: PC_PFMG *ex = (PC_PFMG *)pc->data;
2334: PetscFunctionBegin;
2335: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
2336: if (iascii) {
2337: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE PFMG preconditioning\n"));
2338: PetscCall(PetscViewerASCIIPrintf(viewer, " max iterations %" PetscInt_FMT "\n", ex->its));
2339: PetscCall(PetscViewerASCIIPrintf(viewer, " tolerance %g\n", ex->tol));
2340: PetscCall(PetscViewerASCIIPrintf(viewer, " relax type %s\n", PFMGRelaxType[ex->relax_type]));
2341: PetscCall(PetscViewerASCIIPrintf(viewer, " RAP type %s\n", PFMGRAPType[ex->rap_type]));
2342: PetscCall(PetscViewerASCIIPrintf(viewer, " number pre-relax %" PetscInt_FMT " post-relax %" PetscInt_FMT "\n", ex->num_pre_relax, ex->num_post_relax));
2343: PetscCall(PetscViewerASCIIPrintf(viewer, " max levels %" PetscInt_FMT "\n", ex->max_levels));
2344: PetscCall(PetscViewerASCIIPrintf(viewer, " skip relax %" PetscInt_FMT "\n", ex->skip_relax));
2345: }
2346: PetscFunctionReturn(PETSC_SUCCESS);
2347: }
2349: PetscErrorCode PCSetFromOptions_PFMG(PC pc, PetscOptionItems *PetscOptionsObject)
2350: {
2351: PC_PFMG *ex = (PC_PFMG *)pc->data;
2353: PetscFunctionBegin;
2354: PetscOptionsHeadBegin(PetscOptionsObject, "PFMG options");
2355: PetscCall(PetscOptionsBool("-pc_pfmg_print_statistics", "Print statistics", "HYPRE_StructPFMGSetPrintLevel", ex->print_statistics, &ex->print_statistics, NULL));
2356: PetscCall(PetscOptionsInt("-pc_pfmg_its", "Number of iterations of PFMG to use as preconditioner", "HYPRE_StructPFMGSetMaxIter", ex->its, &ex->its, NULL));
2357: PetscCallExternal(HYPRE_StructPFMGSetMaxIter, ex->hsolver, ex->its);
2358: PetscCall(PetscOptionsInt("-pc_pfmg_num_pre_relax", "Number of smoothing steps before coarse grid", "HYPRE_StructPFMGSetNumPreRelax", ex->num_pre_relax, &ex->num_pre_relax, NULL));
2359: PetscCallExternal(HYPRE_StructPFMGSetNumPreRelax, ex->hsolver, ex->num_pre_relax);
2360: PetscCall(PetscOptionsInt("-pc_pfmg_num_post_relax", "Number of smoothing steps after coarse grid", "HYPRE_StructPFMGSetNumPostRelax", ex->num_post_relax, &ex->num_post_relax, NULL));
2361: PetscCallExternal(HYPRE_StructPFMGSetNumPostRelax, ex->hsolver, ex->num_post_relax);
2363: PetscCall(PetscOptionsInt("-pc_pfmg_max_levels", "Max Levels for MG hierarchy", "HYPRE_StructPFMGSetMaxLevels", ex->max_levels, &ex->max_levels, NULL));
2364: PetscCallExternal(HYPRE_StructPFMGSetMaxLevels, ex->hsolver, ex->max_levels);
2366: PetscCall(PetscOptionsReal("-pc_pfmg_tol", "Tolerance of PFMG", "HYPRE_StructPFMGSetTol", ex->tol, &ex->tol, NULL));
2367: PetscCallExternal(HYPRE_StructPFMGSetTol, ex->hsolver, ex->tol);
2368: PetscCall(PetscOptionsEList("-pc_pfmg_relax_type", "Relax type for the up and down cycles", "HYPRE_StructPFMGSetRelaxType", PFMGRelaxType, PETSC_STATIC_ARRAY_LENGTH(PFMGRelaxType), PFMGRelaxType[ex->relax_type], &ex->relax_type, NULL));
2369: PetscCallExternal(HYPRE_StructPFMGSetRelaxType, ex->hsolver, ex->relax_type);
2370: PetscCall(PetscOptionsEList("-pc_pfmg_rap_type", "RAP type", "HYPRE_StructPFMGSetRAPType", PFMGRAPType, PETSC_STATIC_ARRAY_LENGTH(PFMGRAPType), PFMGRAPType[ex->rap_type], &ex->rap_type, NULL));
2371: PetscCallExternal(HYPRE_StructPFMGSetRAPType, ex->hsolver, ex->rap_type);
2372: PetscCall(PetscOptionsInt("-pc_pfmg_skip_relax", "Skip relaxation on certain grids for isotropic problems. This can greatly improve efficiency by eliminating unnecessary relaxations when the underlying problem is isotropic", "HYPRE_StructPFMGSetSkipRelax", ex->skip_relax, &ex->skip_relax, NULL));
2373: PetscCallExternal(HYPRE_StructPFMGSetSkipRelax, ex->hsolver, ex->skip_relax);
2374: PetscOptionsHeadEnd();
2375: PetscFunctionReturn(PETSC_SUCCESS);
2376: }
2378: PetscErrorCode PCApply_PFMG(PC pc, Vec x, Vec y)
2379: {
2380: PC_PFMG *ex = (PC_PFMG *)pc->data;
2381: PetscScalar *yy;
2382: const PetscScalar *xx;
2383: PetscInt ilower[3], iupper[3];
2384: HYPRE_Int hlower[3], hupper[3];
2385: Mat_HYPREStruct *mx = (Mat_HYPREStruct *)(pc->pmat->data);
2387: PetscFunctionBegin;
2388: PetscCall(PetscCitationsRegister(hypreCitation, &cite));
2389: PetscCall(DMDAGetCorners(mx->da, &ilower[0], &ilower[1], &ilower[2], &iupper[0], &iupper[1], &iupper[2]));
2390: /* when HYPRE_MIXEDINT is defined, sizeof(HYPRE_Int) == 32 */
2391: iupper[0] += ilower[0] - 1;
2392: iupper[1] += ilower[1] - 1;
2393: iupper[2] += ilower[2] - 1;
2394: hlower[0] = (HYPRE_Int)ilower[0];
2395: hlower[1] = (HYPRE_Int)ilower[1];
2396: hlower[2] = (HYPRE_Int)ilower[2];
2397: hupper[0] = (HYPRE_Int)iupper[0];
2398: hupper[1] = (HYPRE_Int)iupper[1];
2399: hupper[2] = (HYPRE_Int)iupper[2];
2401: /* copy x values over to hypre */
2402: PetscCallExternal(HYPRE_StructVectorSetConstantValues, mx->hb, 0.0);
2403: PetscCall(VecGetArrayRead(x, &xx));
2404: PetscCallExternal(HYPRE_StructVectorSetBoxValues, mx->hb, hlower, hupper, (HYPRE_Complex *)xx);
2405: PetscCall(VecRestoreArrayRead(x, &xx));
2406: PetscCallExternal(HYPRE_StructVectorAssemble, mx->hb);
2407: PetscCallExternal(HYPRE_StructPFMGSolve, ex->hsolver, mx->hmat, mx->hb, mx->hx);
2409: /* copy solution values back to PETSc */
2410: PetscCall(VecGetArray(y, &yy));
2411: PetscCallExternal(HYPRE_StructVectorGetBoxValues, mx->hx, hlower, hupper, (HYPRE_Complex *)yy);
2412: PetscCall(VecRestoreArray(y, &yy));
2413: PetscFunctionReturn(PETSC_SUCCESS);
2414: }
2416: 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)
2417: {
2418: PC_PFMG *jac = (PC_PFMG *)pc->data;
2419: HYPRE_Int oits;
2421: PetscFunctionBegin;
2422: PetscCall(PetscCitationsRegister(hypreCitation, &cite));
2423: PetscCallExternal(HYPRE_StructPFMGSetMaxIter, jac->hsolver, its * jac->its);
2424: PetscCallExternal(HYPRE_StructPFMGSetTol, jac->hsolver, rtol);
2426: PetscCall(PCApply_PFMG(pc, b, y));
2427: PetscCallExternal(HYPRE_StructPFMGGetNumIterations, jac->hsolver, &oits);
2428: *outits = oits;
2429: if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
2430: else *reason = PCRICHARDSON_CONVERGED_RTOL;
2431: PetscCallExternal(HYPRE_StructPFMGSetTol, jac->hsolver, jac->tol);
2432: PetscCallExternal(HYPRE_StructPFMGSetMaxIter, jac->hsolver, jac->its);
2433: PetscFunctionReturn(PETSC_SUCCESS);
2434: }
2436: PetscErrorCode PCSetUp_PFMG(PC pc)
2437: {
2438: PC_PFMG *ex = (PC_PFMG *)pc->data;
2439: Mat_HYPREStruct *mx = (Mat_HYPREStruct *)(pc->pmat->data);
2440: PetscBool flg;
2442: PetscFunctionBegin;
2443: PetscCall(PetscObjectTypeCompare((PetscObject)pc->pmat, MATHYPRESTRUCT, &flg));
2444: PetscCheck(flg, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_INCOMP, "Must use MATHYPRESTRUCT with this preconditioner");
2446: /* create the hypre solver object and set its information */
2447: if (ex->hsolver) PetscCallExternal(HYPRE_StructPFMGDestroy, ex->hsolver);
2448: PetscCallExternal(HYPRE_StructPFMGCreate, ex->hcomm, &ex->hsolver);
2450: // Print Hypre statistics about the solve process
2451: if (ex->print_statistics) PetscCallExternal(HYPRE_StructPFMGSetPrintLevel, ex->hsolver, 3);
2453: // The hypre options must be repeated here because the StructPFMG was destroyed and recreated
2454: PetscCallExternal(HYPRE_StructPFMGSetMaxIter, ex->hsolver, ex->its);
2455: PetscCallExternal(HYPRE_StructPFMGSetNumPreRelax, ex->hsolver, ex->num_pre_relax);
2456: PetscCallExternal(HYPRE_StructPFMGSetNumPostRelax, ex->hsolver, ex->num_post_relax);
2457: PetscCallExternal(HYPRE_StructPFMGSetMaxLevels, ex->hsolver, ex->max_levels);
2458: PetscCallExternal(HYPRE_StructPFMGSetTol, ex->hsolver, ex->tol);
2459: PetscCallExternal(HYPRE_StructPFMGSetRelaxType, ex->hsolver, ex->relax_type);
2460: PetscCallExternal(HYPRE_StructPFMGSetRAPType, ex->hsolver, ex->rap_type);
2462: PetscCallExternal(HYPRE_StructPFMGSetup, ex->hsolver, mx->hmat, mx->hb, mx->hx);
2463: PetscCallExternal(HYPRE_StructPFMGSetZeroGuess, ex->hsolver);
2464: PetscFunctionReturn(PETSC_SUCCESS);
2465: }
2467: /*MC
2468: PCPFMG - the hypre PFMG multigrid solver
2470: Options Database Keys:
2471: + -pc_pfmg_its <its> - number of iterations of PFMG to use as preconditioner
2472: . -pc_pfmg_num_pre_relax <steps> - number of smoothing steps before coarse grid solve
2473: . -pc_pfmg_num_post_relax <steps> - number of smoothing steps after coarse grid solve
2474: . -pc_pfmg_tol <tol> - tolerance of PFMG
2475: . -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
2476: . -pc_pfmg_rap_type - type of coarse matrix generation, one of Galerkin,non-Galerkin
2477: - -pc_pfmg_skip_relax - skip relaxation on certain grids for isotropic problems. This can greatly improve efficiency by eliminating unnecessary relaxations
2478: when the underlying problem is isotropic, one of 0,1
2480: Level: advanced
2482: Notes:
2483: This is for CELL-centered descretizations
2485: See `PCSYSPFMG` for a version suitable for systems of PDEs, and `PCSMG`
2487: See `PCHYPRE` for hypre's BoomerAMG algebraic multigrid solver
2489: This must be used with the `MATHYPRESTRUCT` matrix type.
2491: This provides only some of the functionality of PFMG, it supports only one block per process defined by a PETSc `DMDA`.
2493: .seealso: `PCMG`, `MATHYPRESTRUCT`, `PCHYPRE`, `PCGAMG`, `PCSYSPFMG`, `PCSMG`
2494: M*/
2496: PETSC_EXTERN PetscErrorCode PCCreate_PFMG(PC pc)
2497: {
2498: PC_PFMG *ex;
2500: PetscFunctionBegin;
2501: PetscCall(PetscNew(&ex));
2502: pc->data = ex;
2504: ex->its = 1;
2505: ex->tol = 1.e-8;
2506: ex->relax_type = 1;
2507: ex->rap_type = 0;
2508: ex->num_pre_relax = 1;
2509: ex->num_post_relax = 1;
2510: ex->max_levels = 0;
2511: ex->skip_relax = 0;
2512: ex->print_statistics = PETSC_FALSE;
2514: pc->ops->setfromoptions = PCSetFromOptions_PFMG;
2515: pc->ops->view = PCView_PFMG;
2516: pc->ops->destroy = PCDestroy_PFMG;
2517: pc->ops->apply = PCApply_PFMG;
2518: pc->ops->applyrichardson = PCApplyRichardson_PFMG;
2519: pc->ops->setup = PCSetUp_PFMG;
2521: PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
2522: PetscCallExternal(HYPRE_StructPFMGCreate, ex->hcomm, &ex->hsolver);
2523: PetscFunctionReturn(PETSC_SUCCESS);
2524: }
2526: /* we know we are working with a HYPRE_SStructMatrix */
2527: typedef struct {
2528: MPI_Comm hcomm; /* does not share comm with HYPRE_SStructMatrix because need to create solver before getting matrix */
2529: HYPRE_SStructSolver ss_solver;
2531: /* keep copy of SYSPFMG options used so may view them */
2532: PetscInt its;
2533: double tol;
2534: PetscInt relax_type;
2535: PetscInt num_pre_relax, num_post_relax;
2536: } PC_SysPFMG;
2538: PetscErrorCode PCDestroy_SysPFMG(PC pc)
2539: {
2540: PC_SysPFMG *ex = (PC_SysPFMG *)pc->data;
2542: PetscFunctionBegin;
2543: if (ex->ss_solver) PetscCallExternal(HYPRE_SStructSysPFMGDestroy, ex->ss_solver);
2544: PetscCall(PetscCommRestoreComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
2545: PetscCall(PetscFree(pc->data));
2546: PetscFunctionReturn(PETSC_SUCCESS);
2547: }
2549: static const char *SysPFMGRelaxType[] = {"Weighted-Jacobi", "Red/Black-Gauss-Seidel"};
2551: PetscErrorCode PCView_SysPFMG(PC pc, PetscViewer viewer)
2552: {
2553: PetscBool iascii;
2554: PC_SysPFMG *ex = (PC_SysPFMG *)pc->data;
2556: PetscFunctionBegin;
2557: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
2558: if (iascii) {
2559: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE SysPFMG preconditioning\n"));
2560: PetscCall(PetscViewerASCIIPrintf(viewer, " max iterations %" PetscInt_FMT "\n", ex->its));
2561: PetscCall(PetscViewerASCIIPrintf(viewer, " tolerance %g\n", ex->tol));
2562: PetscCall(PetscViewerASCIIPrintf(viewer, " relax type %s\n", PFMGRelaxType[ex->relax_type]));
2563: PetscCall(PetscViewerASCIIPrintf(viewer, " number pre-relax %" PetscInt_FMT " post-relax %" PetscInt_FMT "\n", ex->num_pre_relax, ex->num_post_relax));
2564: }
2565: PetscFunctionReturn(PETSC_SUCCESS);
2566: }
2568: PetscErrorCode PCSetFromOptions_SysPFMG(PC pc, PetscOptionItems *PetscOptionsObject)
2569: {
2570: PC_SysPFMG *ex = (PC_SysPFMG *)pc->data;
2571: PetscBool flg = PETSC_FALSE;
2573: PetscFunctionBegin;
2574: PetscOptionsHeadBegin(PetscOptionsObject, "SysPFMG options");
2575: PetscCall(PetscOptionsBool("-pc_syspfmg_print_statistics", "Print statistics", "HYPRE_SStructSysPFMGSetPrintLevel", flg, &flg, NULL));
2576: if (flg) PetscCallExternal(HYPRE_SStructSysPFMGSetPrintLevel, ex->ss_solver, 3);
2577: PetscCall(PetscOptionsInt("-pc_syspfmg_its", "Number of iterations of SysPFMG to use as preconditioner", "HYPRE_SStructSysPFMGSetMaxIter", ex->its, &ex->its, NULL));
2578: PetscCallExternal(HYPRE_SStructSysPFMGSetMaxIter, ex->ss_solver, ex->its);
2579: PetscCall(PetscOptionsInt("-pc_syspfmg_num_pre_relax", "Number of smoothing steps before coarse grid", "HYPRE_SStructSysPFMGSetNumPreRelax", ex->num_pre_relax, &ex->num_pre_relax, NULL));
2580: PetscCallExternal(HYPRE_SStructSysPFMGSetNumPreRelax, ex->ss_solver, ex->num_pre_relax);
2581: PetscCall(PetscOptionsInt("-pc_syspfmg_num_post_relax", "Number of smoothing steps after coarse grid", "HYPRE_SStructSysPFMGSetNumPostRelax", ex->num_post_relax, &ex->num_post_relax, NULL));
2582: PetscCallExternal(HYPRE_SStructSysPFMGSetNumPostRelax, ex->ss_solver, ex->num_post_relax);
2584: PetscCall(PetscOptionsReal("-pc_syspfmg_tol", "Tolerance of SysPFMG", "HYPRE_SStructSysPFMGSetTol", ex->tol, &ex->tol, NULL));
2585: PetscCallExternal(HYPRE_SStructSysPFMGSetTol, ex->ss_solver, ex->tol);
2586: PetscCall(PetscOptionsEList("-pc_syspfmg_relax_type", "Relax type for the up and down cycles", "HYPRE_SStructSysPFMGSetRelaxType", SysPFMGRelaxType, PETSC_STATIC_ARRAY_LENGTH(SysPFMGRelaxType), SysPFMGRelaxType[ex->relax_type], &ex->relax_type, NULL));
2587: PetscCallExternal(HYPRE_SStructSysPFMGSetRelaxType, ex->ss_solver, ex->relax_type);
2588: PetscOptionsHeadEnd();
2589: PetscFunctionReturn(PETSC_SUCCESS);
2590: }
2592: PetscErrorCode PCApply_SysPFMG(PC pc, Vec x, Vec y)
2593: {
2594: PC_SysPFMG *ex = (PC_SysPFMG *)pc->data;
2595: PetscScalar *yy;
2596: const PetscScalar *xx;
2597: PetscInt ilower[3], iupper[3];
2598: HYPRE_Int hlower[3], hupper[3];
2599: Mat_HYPRESStruct *mx = (Mat_HYPRESStruct *)(pc->pmat->data);
2600: PetscInt ordering = mx->dofs_order;
2601: PetscInt nvars = mx->nvars;
2602: PetscInt part = 0;
2603: PetscInt size;
2604: PetscInt i;
2606: PetscFunctionBegin;
2607: PetscCall(PetscCitationsRegister(hypreCitation, &cite));
2608: PetscCall(DMDAGetCorners(mx->da, &ilower[0], &ilower[1], &ilower[2], &iupper[0], &iupper[1], &iupper[2]));
2609: /* when HYPRE_MIXEDINT is defined, sizeof(HYPRE_Int) == 32 */
2610: iupper[0] += ilower[0] - 1;
2611: iupper[1] += ilower[1] - 1;
2612: iupper[2] += ilower[2] - 1;
2613: hlower[0] = (HYPRE_Int)ilower[0];
2614: hlower[1] = (HYPRE_Int)ilower[1];
2615: hlower[2] = (HYPRE_Int)ilower[2];
2616: hupper[0] = (HYPRE_Int)iupper[0];
2617: hupper[1] = (HYPRE_Int)iupper[1];
2618: hupper[2] = (HYPRE_Int)iupper[2];
2620: size = 1;
2621: for (i = 0; i < 3; i++) size *= (iupper[i] - ilower[i] + 1);
2623: /* copy x values over to hypre for variable ordering */
2624: if (ordering) {
2625: PetscCallExternal(HYPRE_SStructVectorSetConstantValues, mx->ss_b, 0.0);
2626: PetscCall(VecGetArrayRead(x, &xx));
2627: for (i = 0; i < nvars; i++) PetscCallExternal(HYPRE_SStructVectorSetBoxValues, mx->ss_b, part, hlower, hupper, i, (HYPRE_Complex *)(xx + (size * i)));
2628: PetscCall(VecRestoreArrayRead(x, &xx));
2629: PetscCallExternal(HYPRE_SStructVectorAssemble, mx->ss_b);
2630: PetscCallExternal(HYPRE_SStructMatrixMatvec, 1.0, mx->ss_mat, mx->ss_b, 0.0, mx->ss_x);
2631: PetscCallExternal(HYPRE_SStructSysPFMGSolve, ex->ss_solver, mx->ss_mat, mx->ss_b, mx->ss_x);
2633: /* copy solution values back to PETSc */
2634: PetscCall(VecGetArray(y, &yy));
2635: for (i = 0; i < nvars; i++) PetscCallExternal(HYPRE_SStructVectorGetBoxValues, mx->ss_x, part, hlower, hupper, i, (HYPRE_Complex *)(yy + (size * i)));
2636: PetscCall(VecRestoreArray(y, &yy));
2637: } else { /* nodal ordering must be mapped to variable ordering for sys_pfmg */
2638: PetscScalar *z;
2639: PetscInt j, k;
2641: PetscCall(PetscMalloc1(nvars * size, &z));
2642: PetscCallExternal(HYPRE_SStructVectorSetConstantValues, mx->ss_b, 0.0);
2643: PetscCall(VecGetArrayRead(x, &xx));
2645: /* transform nodal to hypre's variable ordering for sys_pfmg */
2646: for (i = 0; i < size; i++) {
2647: k = i * nvars;
2648: for (j = 0; j < nvars; j++) z[j * size + i] = xx[k + j];
2649: }
2650: for (i = 0; i < nvars; i++) PetscCallExternal(HYPRE_SStructVectorSetBoxValues, mx->ss_b, part, hlower, hupper, i, (HYPRE_Complex *)(z + (size * i)));
2651: PetscCall(VecRestoreArrayRead(x, &xx));
2652: PetscCallExternal(HYPRE_SStructVectorAssemble, mx->ss_b);
2653: PetscCallExternal(HYPRE_SStructSysPFMGSolve, ex->ss_solver, mx->ss_mat, mx->ss_b, mx->ss_x);
2655: /* copy solution values back to PETSc */
2656: PetscCall(VecGetArray(y, &yy));
2657: for (i = 0; i < nvars; i++) PetscCallExternal(HYPRE_SStructVectorGetBoxValues, mx->ss_x, part, hlower, hupper, i, (HYPRE_Complex *)(z + (size * i)));
2658: /* transform hypre's variable ordering for sys_pfmg to nodal ordering */
2659: for (i = 0; i < size; i++) {
2660: k = i * nvars;
2661: for (j = 0; j < nvars; j++) yy[k + j] = z[j * size + i];
2662: }
2663: PetscCall(VecRestoreArray(y, &yy));
2664: PetscCall(PetscFree(z));
2665: }
2666: PetscFunctionReturn(PETSC_SUCCESS);
2667: }
2669: 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)
2670: {
2671: PC_SysPFMG *jac = (PC_SysPFMG *)pc->data;
2672: HYPRE_Int oits;
2674: PetscFunctionBegin;
2675: PetscCall(PetscCitationsRegister(hypreCitation, &cite));
2676: PetscCallExternal(HYPRE_SStructSysPFMGSetMaxIter, jac->ss_solver, its * jac->its);
2677: PetscCallExternal(HYPRE_SStructSysPFMGSetTol, jac->ss_solver, rtol);
2678: PetscCall(PCApply_SysPFMG(pc, b, y));
2679: PetscCallExternal(HYPRE_SStructSysPFMGGetNumIterations, jac->ss_solver, &oits);
2680: *outits = oits;
2681: if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
2682: else *reason = PCRICHARDSON_CONVERGED_RTOL;
2683: PetscCallExternal(HYPRE_SStructSysPFMGSetTol, jac->ss_solver, jac->tol);
2684: PetscCallExternal(HYPRE_SStructSysPFMGSetMaxIter, jac->ss_solver, jac->its);
2685: PetscFunctionReturn(PETSC_SUCCESS);
2686: }
2688: PetscErrorCode PCSetUp_SysPFMG(PC pc)
2689: {
2690: PC_SysPFMG *ex = (PC_SysPFMG *)pc->data;
2691: Mat_HYPRESStruct *mx = (Mat_HYPRESStruct *)(pc->pmat->data);
2692: PetscBool flg;
2694: PetscFunctionBegin;
2695: PetscCall(PetscObjectTypeCompare((PetscObject)pc->pmat, MATHYPRESSTRUCT, &flg));
2696: PetscCheck(flg, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_INCOMP, "Must use MATHYPRESSTRUCT with this preconditioner");
2698: /* create the hypre sstruct solver object and set its information */
2699: if (ex->ss_solver) PetscCallExternal(HYPRE_SStructSysPFMGDestroy, ex->ss_solver);
2700: PetscCallExternal(HYPRE_SStructSysPFMGCreate, ex->hcomm, &ex->ss_solver);
2701: PetscCallExternal(HYPRE_SStructSysPFMGSetZeroGuess, ex->ss_solver);
2702: PetscCallExternal(HYPRE_SStructSysPFMGSetup, ex->ss_solver, mx->ss_mat, mx->ss_b, mx->ss_x);
2703: PetscFunctionReturn(PETSC_SUCCESS);
2704: }
2706: /*MC
2707: PCSYSPFMG - the hypre SysPFMG multigrid solver
2709: Level: advanced
2711: Options Database Keys:
2712: + -pc_syspfmg_its <its> - number of iterations of SysPFMG to use as preconditioner
2713: . -pc_syspfmg_num_pre_relax <steps> - number of smoothing steps before coarse grid
2714: . -pc_syspfmg_num_post_relax <steps> - number of smoothing steps after coarse grid
2715: . -pc_syspfmg_tol <tol> - tolerance of SysPFMG
2716: - -pc_syspfmg_relax_type <Weighted-Jacobi,Red/Black-Gauss-Seidel> - relaxation type for the up and down cycles
2718: Notes:
2719: See `PCPFMG` for hypre's PFMG that works for a scalar PDE and `PCSMG`
2721: See `PCHYPRE` for hypre's BoomerAMG algebraic multigrid solver
2723: This is for CELL-centered descretizations
2725: This must be used with the `MATHYPRESSTRUCT` matrix type.
2727: This does not give access to all the functionality of hypres SysPFMG, it supports only one part, and one block per process defined by a PETSc `DMDA`.
2729: .seealso: `PCMG`, `MATHYPRESSTRUCT`, `PCPFMG`, `PCHYPRE`, `PCGAMG`, `PCSMG`
2730: M*/
2732: PETSC_EXTERN PetscErrorCode PCCreate_SysPFMG(PC pc)
2733: {
2734: PC_SysPFMG *ex;
2736: PetscFunctionBegin;
2737: PetscCall(PetscNew(&ex));
2738: pc->data = ex;
2740: ex->its = 1;
2741: ex->tol = 1.e-8;
2742: ex->relax_type = 1;
2743: ex->num_pre_relax = 1;
2744: ex->num_post_relax = 1;
2746: pc->ops->setfromoptions = PCSetFromOptions_SysPFMG;
2747: pc->ops->view = PCView_SysPFMG;
2748: pc->ops->destroy = PCDestroy_SysPFMG;
2749: pc->ops->apply = PCApply_SysPFMG;
2750: pc->ops->applyrichardson = PCApplyRichardson_SysPFMG;
2751: pc->ops->setup = PCSetUp_SysPFMG;
2753: PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
2754: PetscCallExternal(HYPRE_SStructSysPFMGCreate, ex->hcomm, &ex->ss_solver);
2755: PetscFunctionReturn(PETSC_SUCCESS);
2756: }
2758: /* PC SMG */
2759: typedef struct {
2760: MPI_Comm hcomm; /* does not share comm with HYPRE_StructMatrix because need to create solver before getting matrix */
2761: HYPRE_StructSolver hsolver;
2762: PetscInt its; /* keep copy of SMG options used so may view them */
2763: double tol;
2764: PetscBool print_statistics;
2765: PetscInt num_pre_relax, num_post_relax;
2766: } PC_SMG;
2768: PetscErrorCode PCDestroy_SMG(PC pc)
2769: {
2770: PC_SMG *ex = (PC_SMG *)pc->data;
2772: PetscFunctionBegin;
2773: if (ex->hsolver) PetscCallExternal(HYPRE_StructSMGDestroy, ex->hsolver);
2774: PetscCall(PetscCommRestoreComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
2775: PetscCall(PetscFree(pc->data));
2776: PetscFunctionReturn(PETSC_SUCCESS);
2777: }
2779: PetscErrorCode PCView_SMG(PC pc, PetscViewer viewer)
2780: {
2781: PetscBool iascii;
2782: PC_SMG *ex = (PC_SMG *)pc->data;
2784: PetscFunctionBegin;
2785: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
2786: if (iascii) {
2787: PetscCall(PetscViewerASCIIPrintf(viewer, " HYPRE SMG preconditioning\n"));
2788: PetscCall(PetscViewerASCIIPrintf(viewer, " max iterations %" PetscInt_FMT "\n", ex->its));
2789: PetscCall(PetscViewerASCIIPrintf(viewer, " tolerance %g\n", ex->tol));
2790: PetscCall(PetscViewerASCIIPrintf(viewer, " number pre-relax %" PetscInt_FMT " post-relax %" PetscInt_FMT "\n", ex->num_pre_relax, ex->num_post_relax));
2791: }
2792: PetscFunctionReturn(PETSC_SUCCESS);
2793: }
2795: PetscErrorCode PCSetFromOptions_SMG(PC pc, PetscOptionItems *PetscOptionsObject)
2796: {
2797: PC_SMG *ex = (PC_SMG *)pc->data;
2799: PetscFunctionBegin;
2800: PetscOptionsHeadBegin(PetscOptionsObject, "SMG options");
2802: PetscCall(PetscOptionsInt("-pc_smg_its", "Number of iterations of SMG to use as preconditioner", "HYPRE_StructSMGSetMaxIter", ex->its, &ex->its, NULL));
2803: PetscCall(PetscOptionsInt("-pc_smg_num_pre_relax", "Number of smoothing steps before coarse grid", "HYPRE_StructSMGSetNumPreRelax", ex->num_pre_relax, &ex->num_pre_relax, NULL));
2804: PetscCall(PetscOptionsInt("-pc_smg_num_post_relax", "Number of smoothing steps after coarse grid", "HYPRE_StructSMGSetNumPostRelax", ex->num_post_relax, &ex->num_post_relax, NULL));
2805: PetscCall(PetscOptionsReal("-pc_smg_tol", "Tolerance of SMG", "HYPRE_StructSMGSetTol", ex->tol, &ex->tol, NULL));
2807: PetscOptionsHeadEnd();
2808: PetscFunctionReturn(PETSC_SUCCESS);
2809: }
2811: PetscErrorCode PCApply_SMG(PC pc, Vec x, Vec y)
2812: {
2813: PC_SMG *ex = (PC_SMG *)pc->data;
2814: PetscScalar *yy;
2815: const PetscScalar *xx;
2816: PetscInt ilower[3], iupper[3];
2817: HYPRE_Int hlower[3], hupper[3];
2818: Mat_HYPREStruct *mx = (Mat_HYPREStruct *)(pc->pmat->data);
2820: PetscFunctionBegin;
2821: PetscCall(PetscCitationsRegister(hypreCitation, &cite));
2822: PetscCall(DMDAGetCorners(mx->da, &ilower[0], &ilower[1], &ilower[2], &iupper[0], &iupper[1], &iupper[2]));
2823: /* when HYPRE_MIXEDINT is defined, sizeof(HYPRE_Int) == 32 */
2824: iupper[0] += ilower[0] - 1;
2825: iupper[1] += ilower[1] - 1;
2826: iupper[2] += ilower[2] - 1;
2827: hlower[0] = (HYPRE_Int)ilower[0];
2828: hlower[1] = (HYPRE_Int)ilower[1];
2829: hlower[2] = (HYPRE_Int)ilower[2];
2830: hupper[0] = (HYPRE_Int)iupper[0];
2831: hupper[1] = (HYPRE_Int)iupper[1];
2832: hupper[2] = (HYPRE_Int)iupper[2];
2834: /* copy x values over to hypre */
2835: PetscCallExternal(HYPRE_StructVectorSetConstantValues, mx->hb, 0.0);
2836: PetscCall(VecGetArrayRead(x, &xx));
2837: PetscCallExternal(HYPRE_StructVectorSetBoxValues, mx->hb, hlower, hupper, (HYPRE_Complex *)xx);
2838: PetscCall(VecRestoreArrayRead(x, &xx));
2839: PetscCallExternal(HYPRE_StructVectorAssemble, mx->hb);
2840: PetscCallExternal(HYPRE_StructSMGSolve, ex->hsolver, mx->hmat, mx->hb, mx->hx);
2842: /* copy solution values back to PETSc */
2843: PetscCall(VecGetArray(y, &yy));
2844: PetscCallExternal(HYPRE_StructVectorGetBoxValues, mx->hx, hlower, hupper, (HYPRE_Complex *)yy);
2845: PetscCall(VecRestoreArray(y, &yy));
2846: PetscFunctionReturn(PETSC_SUCCESS);
2847: }
2849: static PetscErrorCode PCApplyRichardson_SMG(PC pc, Vec b, Vec y, Vec w, PetscReal rtol, PetscReal abstol, PetscReal dtol, PetscInt its, PetscBool guesszero, PetscInt *outits, PCRichardsonConvergedReason *reason)
2850: {
2851: PC_SMG *jac = (PC_SMG *)pc->data;
2852: HYPRE_Int oits;
2854: PetscFunctionBegin;
2855: PetscCall(PetscCitationsRegister(hypreCitation, &cite));
2856: PetscCallExternal(HYPRE_StructSMGSetMaxIter, jac->hsolver, its * jac->its);
2857: PetscCallExternal(HYPRE_StructSMGSetTol, jac->hsolver, rtol);
2859: PetscCall(PCApply_SMG(pc, b, y));
2860: PetscCallExternal(HYPRE_StructSMGGetNumIterations, jac->hsolver, &oits);
2861: *outits = oits;
2862: if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
2863: else *reason = PCRICHARDSON_CONVERGED_RTOL;
2864: PetscCallExternal(HYPRE_StructSMGSetTol, jac->hsolver, jac->tol);
2865: PetscCallExternal(HYPRE_StructSMGSetMaxIter, jac->hsolver, jac->its);
2866: PetscFunctionReturn(PETSC_SUCCESS);
2867: }
2869: PetscErrorCode PCSetUp_SMG(PC pc)
2870: {
2871: PetscInt i, dim;
2872: PC_SMG *ex = (PC_SMG *)pc->data;
2873: Mat_HYPREStruct *mx = (Mat_HYPREStruct *)(pc->pmat->data);
2874: PetscBool flg;
2875: DMBoundaryType p[3];
2876: PetscInt M[3];
2878: PetscFunctionBegin;
2879: PetscCall(PetscObjectTypeCompare((PetscObject)pc->pmat, MATHYPRESTRUCT, &flg));
2880: PetscCheck(flg, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_INCOMP, "Must use MATHYPRESTRUCT with this preconditioner");
2882: PetscCall(DMDAGetInfo(mx->da, &dim, &M[0], &M[1], &M[2], 0, 0, 0, 0, 0, &p[0], &p[1], &p[2], 0));
2883: // Check if power of 2 in periodic directions
2884: for (i = 0; i < dim; i++) {
2885: if (((M[i] & (M[i] - 1)) != 0) && (p[i] == DM_BOUNDARY_PERIODIC)) {
2886: SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_INCOMP, "With SMG, the number of points in a periodic direction must be a power of 2, but is here %" PetscInt_FMT ".", M[i]);
2887: }
2888: }
2890: /* create the hypre solver object and set its information */
2891: if (ex->hsolver) PetscCallExternal(HYPRE_StructSMGDestroy, (ex->hsolver));
2892: PetscCallExternal(HYPRE_StructSMGCreate, ex->hcomm, &ex->hsolver);
2893: // The hypre options must be set here and not in SetFromOptions because it is created here!
2894: PetscCallExternal(HYPRE_StructSMGSetMaxIter, ex->hsolver, ex->its);
2895: PetscCallExternal(HYPRE_StructSMGSetNumPreRelax, ex->hsolver, ex->num_pre_relax);
2896: PetscCallExternal(HYPRE_StructSMGSetNumPostRelax, ex->hsolver, ex->num_post_relax);
2897: PetscCallExternal(HYPRE_StructSMGSetTol, ex->hsolver, ex->tol);
2899: PetscCallExternal(HYPRE_StructSMGSetup, ex->hsolver, mx->hmat, mx->hb, mx->hx);
2900: PetscCallExternal(HYPRE_StructSMGSetZeroGuess, ex->hsolver);
2901: PetscFunctionReturn(PETSC_SUCCESS);
2902: }
2904: /*MC
2905: PCSMG - the hypre (structured grid) SMG multigrid solver
2907: Level: advanced
2909: Options Database Keys:
2910: + -pc_smg_its <its> - number of iterations of SMG to use as preconditioner
2911: . -pc_smg_num_pre_relax <steps> - number of smoothing steps before coarse grid
2912: . -pc_smg_num_post_relax <steps> - number of smoothing steps after coarse grid
2913: - -pc_smg_tol <tol> - tolerance of SMG
2915: Notes:
2916: This is for CELL-centered descretizations
2918: This must be used with the `MATHYPRESTRUCT` `MatType`.
2920: This does not provide all the functionality of hypre's SMG solver, it supports only one block per process defined by a PETSc `DMDA`.
2922: See `PCSYSPFMG`, `PCSMG`, `PCPFMG`, and `PCHYPRE` for access to hypre's other preconditioners
2924: .seealso: `PCMG`, `MATHYPRESTRUCT`, `PCPFMG`, `PCSYSPFMG`, `PCHYPRE`, `PCGAMG`
2925: M*/
2927: PETSC_EXTERN PetscErrorCode PCCreate_SMG(PC pc)
2928: {
2929: PC_SMG *ex;
2931: PetscFunctionBegin;
2932: PetscCall(PetscNew(&ex));
2933: pc->data = ex;
2935: ex->its = 1;
2936: ex->tol = 1.e-8;
2937: ex->num_pre_relax = 1;
2938: ex->num_post_relax = 1;
2940: pc->ops->setfromoptions = PCSetFromOptions_SMG;
2941: pc->ops->view = PCView_SMG;
2942: pc->ops->destroy = PCDestroy_SMG;
2943: pc->ops->apply = PCApply_SMG;
2944: pc->ops->applyrichardson = PCApplyRichardson_SMG;
2945: pc->ops->setup = PCSetUp_SMG;
2947: PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
2948: PetscCallExternal(HYPRE_StructSMGCreate, ex->hcomm, &ex->hsolver);
2949: PetscFunctionReturn(PETSC_SUCCESS);
2950: }