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:   PetscObjectParameterDeclare(PetscInt, coarsentype);
 64:   PetscInt  measuretype;
 65:   PetscInt  smoothtype;
 66:   PetscInt  smoothsweeps;
 67:   PetscInt  smoothnumlevels;
 68:   PetscInt  eu_level;         /* Number of levels for ILU(k) in Euclid */
 69:   PetscReal eu_droptolerance; /* Drop tolerance for ILU(k) in Euclid */
 70:   PetscInt  eu_bj;            /* Defines use of Block Jacobi ILU in Euclid */
 71:   PetscObjectParameterDeclare(PetscInt, relaxtype[3]);
 72:   PetscReal relaxweight;
 73:   PetscReal outerrelaxweight;
 74:   PetscObjectParameterDeclare(PetscInt, relaxorder);
 75:   PetscReal truncfactor;
 76:   PetscBool applyrichardson;
 77:   PetscInt  pmax;
 78:   PetscObjectParameterDeclare(PetscInt, interptype);
 79:   PetscInt maxc;
 80:   PetscInt minc;
 81: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
 82:   PetscObjectParameterDeclarePtr(const char, spgemm_type); // this is a global hypre parameter but is closely associated with BoomerAMG
 83: #endif
 84:   /* GPU */
 85:   PetscObjectParameterDeclare(PetscBool, keeptranspose);
 86:   PetscInt rap2;
 87:   PetscObjectParameterDeclare(PetscInt, mod_rap2);

 89:   /* AIR */
 90:   PetscInt  Rtype;
 91:   PetscReal Rstrongthreshold;
 92:   PetscReal Rfilterthreshold;
 93:   PetscInt  Adroptype;
 94:   PetscReal Adroptol;

 96:   PetscInt agg_nl;
 97:   PetscObjectParameterDeclare(PetscInt, agg_interptype);
 98:   PetscInt  agg_num_paths;
 99:   PetscBool nodal_relax;
100:   PetscInt  nodal_relax_levels;

102:   PetscInt  nodal_coarsening;
103:   PetscInt  nodal_coarsening_diag;
104:   PetscInt  vec_interp_variant;
105:   PetscInt  vec_interp_qmax;
106:   PetscBool vec_interp_smooth;
107:   PetscInt  interp_refine;

109:   /* NearNullSpace support */
110:   VecHYPRE_IJVector *hmnull;
111:   HYPRE_ParVector   *phmnull;
112:   PetscInt           n_hmnull;
113:   Vec                hmnull_constant;

115:   /* options for AS (Auxiliary Space preconditioners) */
116:   PetscInt  as_print;
117:   PetscInt  as_max_iter;
118:   PetscReal as_tol;
119:   PetscInt  as_relax_type;
120:   PetscInt  as_relax_times;
121:   PetscReal as_relax_weight;
122:   PetscReal as_omega;
123:   PetscInt  as_amg_alpha_opts[5]; /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for vector Poisson (AMS) or Curl problem (ADS) */
124:   PetscReal as_amg_alpha_theta;   /* AMG strength for vector Poisson (AMS) or Curl problem (ADS) */
125:   PetscInt  as_amg_beta_opts[5];  /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for scalar Poisson (AMS) or vector Poisson (ADS) */
126:   PetscReal as_amg_beta_theta;    /* AMG strength for scalar Poisson (AMS) or vector Poisson (ADS)  */
127:   PetscInt  ams_cycle_type;
128:   PetscInt  ads_cycle_type;

130:   /* additional data */
131:   Mat G;             /* MatHYPRE */
132:   Mat C;             /* MatHYPRE */
133:   Mat alpha_Poisson; /* MatHYPRE */
134:   Mat beta_Poisson;  /* MatHYPRE */

136:   /* extra information for AMS */
137:   PetscInt          dim; /* geometrical dimension */
138:   VecHYPRE_IJVector coords[3];
139:   VecHYPRE_IJVector constants[3];
140:   VecHYPRE_IJVector interior;
141:   Mat               RT_PiFull, RT_Pi[3];
142:   Mat               ND_PiFull, ND_Pi[3];
143:   PetscBool         ams_beta_is_zero;
144:   PetscBool         ams_beta_is_zero_part;
145:   PetscInt          ams_proj_freq;
146: } PC_HYPRE;

148: /*
149:   Matrices with AIJ format are created IN PLACE with using (I,J,data) from BoomerAMG. Since the data format in hypre_ParCSRMatrix
150:   is different from that used in PETSc, the original hypre_ParCSRMatrix can not be used any more after call this routine.
151:   It is used in PCHMG. Other users should avoid using this function.
152: */
153: static PetscErrorCode PCGetCoarseOperators_BoomerAMG(PC pc, PetscInt *nlevels, Mat *operators[])
154: {
155:   PC_HYPRE            *jac = (PC_HYPRE *)pc->data;
156:   PetscBool            same;
157:   PetscInt             num_levels, l;
158:   Mat                 *mattmp;
159:   hypre_ParCSRMatrix **A_array;

161:   PetscFunctionBegin;
162:   PetscCall(PetscStrcmp(jac->hypre_type, "boomeramg", &same));
163:   PetscCheck(same, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_NOTSAMETYPE, "Hypre type is not BoomerAMG");
164:   num_levels = hypre_ParAMGDataNumLevels((hypre_ParAMGData *)jac->hsolver);
165:   PetscCall(PetscMalloc1(num_levels, &mattmp));
166:   A_array = hypre_ParAMGDataAArray((hypre_ParAMGData *)jac->hsolver);
167:   for (l = 1; l < num_levels; l++) {
168:     PetscCall(MatCreateFromParCSR(A_array[l], MATAIJ, PETSC_OWN_POINTER, &mattmp[num_levels - 1 - l]));
169:     /* We want to own the data, and HYPRE can not touch this matrix any more */
170:     A_array[l] = NULL;
171:   }
172:   *nlevels   = num_levels;
173:   *operators = mattmp;
174:   PetscFunctionReturn(PETSC_SUCCESS);
175: }

177: /*
178:   Matrices with AIJ format are created IN PLACE with using (I,J,data) from BoomerAMG. Since the data format in hypre_ParCSRMatrix
179:   is different from that used in PETSc, the original hypre_ParCSRMatrix can not be used any more after call this routine.
180:   It is used in PCHMG. Other users should avoid using this function.
181: */
182: static PetscErrorCode PCGetInterpolations_BoomerAMG(PC pc, PetscInt *nlevels, Mat *interpolations[])
183: {
184:   PC_HYPRE            *jac = (PC_HYPRE *)pc->data;
185:   PetscBool            same;
186:   PetscInt             num_levels, l;
187:   Mat                 *mattmp;
188:   hypre_ParCSRMatrix **P_array;

190:   PetscFunctionBegin;
191:   PetscCall(PetscStrcmp(jac->hypre_type, "boomeramg", &same));
192:   PetscCheck(same, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_NOTSAMETYPE, "Hypre type is not BoomerAMG");
193:   num_levels = hypre_ParAMGDataNumLevels((hypre_ParAMGData *)jac->hsolver);
194:   PetscCall(PetscMalloc1(num_levels, &mattmp));
195:   P_array = hypre_ParAMGDataPArray((hypre_ParAMGData *)jac->hsolver);
196:   for (l = 1; l < num_levels; l++) {
197:     PetscCall(MatCreateFromParCSR(P_array[num_levels - 1 - l], MATAIJ, PETSC_OWN_POINTER, &mattmp[l - 1]));
198:     /* We want to own the data, and HYPRE can not touch this matrix any more */
199:     P_array[num_levels - 1 - l] = NULL;
200:   }
201:   *nlevels        = num_levels;
202:   *interpolations = mattmp;
203:   PetscFunctionReturn(PETSC_SUCCESS);
204: }

206: /*
207:   Boolean Vecs are created IN PLACE with using data from BoomerAMG.
208: */
209: static PetscErrorCode PCHYPREGetCFMarkers_BoomerAMG(PC pc, PetscInt *n_per_level[], PetscBT *CFMarkers[])
210: {
211:   PC_HYPRE        *jac = (PC_HYPRE *)pc->data;
212:   PetscBool        same;
213:   PetscInt         num_levels, fine_nodes = 0, coarse_nodes;
214:   PetscInt        *n_per_temp;
215:   PetscBT         *markertmp;
216:   hypre_IntArray **CF_marker_array;

218:   PetscFunctionBegin;
219:   PetscCall(PetscStrcmp(jac->hypre_type, "boomeramg", &same));
220:   PetscCheck(same, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_NOTSAMETYPE, "Hypre type is not BoomerAMG");
221:   num_levels = hypre_ParAMGDataNumLevels((hypre_ParAMGData *)jac->hsolver);
222:   PetscCall(PetscMalloc1(num_levels, &n_per_temp));
223:   PetscCall(PetscMalloc1(num_levels - 1, &markertmp));
224:   CF_marker_array = hypre_ParAMGDataCFMarkerArray((hypre_ParAMGData *)jac->hsolver);
225:   for (PetscInt l = 0, CFMaxIndex = num_levels - 2; CFMaxIndex >= 0; l++, CFMaxIndex--) {
226:     fine_nodes   = hypre_IntArraySize(CF_marker_array[CFMaxIndex]);
227:     coarse_nodes = 0;
228:     PetscCall(PetscBTCreate(fine_nodes, &markertmp[l]));
229:     for (PetscInt k = 0; k < fine_nodes; k++) {
230:       if (hypre_IntArrayDataI(CF_marker_array[CFMaxIndex], k) > 0) {
231:         PetscCall(PetscBTSet(markertmp[l], k));
232:         coarse_nodes++;
233:       }
234:     }
235:     n_per_temp[l] = coarse_nodes;
236:   }
237:   n_per_temp[num_levels - 1] = fine_nodes;
238:   *n_per_level               = n_per_temp;
239:   *CFMarkers                 = markertmp;
240:   PetscFunctionReturn(PETSC_SUCCESS);
241: }

243: /* Resets (frees) Hypre's representation of the near null space */
244: static PetscErrorCode PCHYPREResetNearNullSpace_Private(PC pc)
245: {
246:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
247:   PetscInt  i;

249:   PetscFunctionBegin;
250:   for (i = 0; i < jac->n_hmnull; i++) PetscCall(VecHYPRE_IJVectorDestroy(&jac->hmnull[i]));
251:   PetscCall(PetscFree(jac->hmnull));
252:   PetscCall(PetscFree(jac->phmnull));
253:   PetscCall(VecDestroy(&jac->hmnull_constant));
254:   jac->n_hmnull = 0;
255:   PetscFunctionReturn(PETSC_SUCCESS);
256: }

258: static const char    *HYPRESpgemmTypes[] = {"cusparse", "hypre"};
259: static PetscErrorCode PCMGGalerkinSetMatProductAlgorithm_HYPRE_BoomerAMG(PC pc, const char name[])
260: {
261:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
262:   PetscBool flag;

264: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
265:   PetscFunctionBegin;
266:   if (jac->spgemm_type) {
267:     PetscCall(PetscStrcmp(jac->spgemm_type, name, &flag));
268:     PetscCheck(flag, PetscObjectComm((PetscObject)pc), PETSC_ERR_ORDER, "PETSc support for resetting the HYPRE SpGEMM is not implemented");
269:     PetscFunctionReturn(PETSC_SUCCESS);
270:   } else jac->spgemm_type = name;

272:   PetscCall(PetscStrcmp("cusparse", jac->spgemm_type, &flag));
273:   if (flag) {
274:     PetscCallExternal(HYPRE_SetSpGemmUseCusparse, 1);
275:     PetscFunctionReturn(PETSC_SUCCESS);
276:   }
277:   PetscCall(PetscStrcmp("hypre", jac->spgemm_type, &flag));
278:   if (flag) {
279:     PetscCallExternal(HYPRE_SetSpGemmUseCusparse, 0);
280:     PetscFunctionReturn(PETSC_SUCCESS);
281:   }
282:   jac->spgemm_type = NULL;
283:   SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown HYPRE SpGEMM type %s; Choices are cusparse, hypre", name);
284: #endif
285: }

287: static PetscErrorCode PCSetUp_HYPRE(PC pc)
288: {
289:   PC_HYPRE          *jac = (PC_HYPRE *)pc->data;
290:   Mat_HYPRE         *hjac;
291:   HYPRE_ParCSRMatrix hmat;
292:   HYPRE_ParVector    bv, xv;
293:   PetscBool          ishypre;

295:   PetscFunctionBegin;
296:   /* default type is boomerAMG */
297:   if (!jac->hypre_type) PetscCall(PCHYPRESetType(pc, "boomeramg"));

299:   /* get hypre matrix */
300:   if (pc->flag == DIFFERENT_NONZERO_PATTERN) PetscCall(MatDestroy(&jac->hpmat));
301:   PetscCall(PetscObjectTypeCompare((PetscObject)pc->pmat, MATHYPRE, &ishypre));
302:   if (!ishypre) {
303: #if defined(PETSC_HAVE_HYPRE_DEVICE) && PETSC_PKG_HYPRE_VERSION_LE(2, 30, 0)
304:     /* Temporary fix since we do not support MAT_REUSE_MATRIX with HYPRE device */
305:     PetscBool iscuda, iship, iskokkos;

307:     PetscCall(PetscObjectTypeCompareAny((PetscObject)pc->pmat, &iscuda, MATSEQAIJCUSPARSE, MATMPIAIJCUSPARSE, ""));
308:     PetscCall(PetscObjectTypeCompareAny((PetscObject)pc->pmat, &iship, MATSEQAIJHIPSPARSE, MATMPIAIJHIPSPARSE, ""));
309:     PetscCall(PetscObjectTypeCompareAny((PetscObject)pc->pmat, &iskokkos, MATSEQAIJKOKKOS, MATMPIAIJKOKKOS, ""));
310:     if (iscuda || iship || iskokkos) PetscCall(MatDestroy(&jac->hpmat));
311: #endif
312:     PetscCall(MatConvert(pc->pmat, MATHYPRE, jac->hpmat ? MAT_REUSE_MATRIX : MAT_INITIAL_MATRIX, &jac->hpmat));
313:   } else {
314:     PetscCall(PetscObjectReference((PetscObject)pc->pmat));
315:     PetscCall(MatDestroy(&jac->hpmat));
316:     jac->hpmat = pc->pmat;
317:   }

319:   /* allow debug */
320:   PetscCall(MatViewFromOptions(jac->hpmat, NULL, "-pc_hypre_mat_view"));
321:   hjac = (Mat_HYPRE *)jac->hpmat->data;

323:   /* special case for BoomerAMG */
324:   if (jac->setup == HYPRE_BoomerAMGSetup) {
325:     MatNullSpace mnull;
326:     PetscBool    has_const;
327:     PetscInt     bs, nvec, i;
328:     PetscMemType memtype;
329:     const Vec   *vecs;

331:     PetscCall(MatGetCurrentMemType(jac->hpmat, &memtype));
332:     if (PetscMemTypeDevice(memtype)) {
333:       /* GPU defaults
334:          From https://hypre.readthedocs.io/en/latest/solvers-boomeramg.html#gpu-supported-options
335:          and /src/parcsr_ls/par_amg.c
336:          First handle options which users have interfaces for changing */
337:       PetscObjectParameterSetDefault(jac, coarsentype, 8);
338:       PetscObjectParameterSetDefault(jac, relaxorder, 0);
339:       PetscObjectParameterSetDefault(jac, interptype, 6);
340:       PetscObjectParameterSetDefault(jac, relaxtype[0], 18);
341:       PetscObjectParameterSetDefault(jac, relaxtype[1], 18);
342: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
343:       PetscObjectParameterSetDefault(jac, spgemm_type, HYPRESpgemmTypes[0]);
344: #endif
345: #if PETSC_PKG_HYPRE_VERSION_GE(2, 18, 0)
346:       PetscObjectParameterSetDefault(jac, keeptranspose, PETSC_TRUE);
347:       PetscObjectParameterSetDefault(jac, mod_rap2, 1);
348: #endif
349:       PetscObjectParameterSetDefault(jac, agg_interptype, 7);
350:     } else {
351:       PetscObjectParameterSetDefault(jac, coarsentype, 6);
352:       PetscObjectParameterSetDefault(jac, relaxorder, 1);
353:       PetscObjectParameterSetDefault(jac, interptype, 0);
354:       PetscObjectParameterSetDefault(jac, relaxtype[0], 6);
355:       PetscObjectParameterSetDefault(jac, relaxtype[1], 6); /* Defaults to SYMMETRIC since in PETSc we are using a PC - most likely with CG */
356: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
357:       PetscObjectParameterSetDefault(jac, spgemm_type, "hypre");
358: #endif
359: #if PETSC_PKG_HYPRE_VERSION_GE(2, 18, 0)
360:       PetscObjectParameterSetDefault(jac, keeptranspose, PETSC_FALSE);
361:       PetscObjectParameterSetDefault(jac, mod_rap2, 0);
362: #endif
363:       PetscObjectParameterSetDefault(jac, agg_interptype, 4);
364:     }
365:     PetscObjectParameterSetDefault(jac, relaxtype[2], 9); /*G.E. */

367:     PetscCallExternal(HYPRE_BoomerAMGSetCycleType, jac->hsolver, jac->cycletype);
368:     PetscCallExternal(HYPRE_BoomerAMGSetMaxLevels, jac->hsolver, jac->maxlevels);
369:     PetscCallExternal(HYPRE_BoomerAMGSetMaxIter, jac->hsolver, jac->maxiter);
370:     PetscCallExternal(HYPRE_BoomerAMGSetTol, jac->hsolver, jac->tol);
371:     PetscCallExternal(HYPRE_BoomerAMGSetTruncFactor, jac->hsolver, jac->truncfactor);
372:     PetscCallExternal(HYPRE_BoomerAMGSetStrongThreshold, jac->hsolver, jac->strongthreshold);
373:     PetscCallExternal(HYPRE_BoomerAMGSetMaxRowSum, jac->hsolver, jac->maxrowsum);
374:     PetscCallExternal(HYPRE_BoomerAMGSetMeasureType, jac->hsolver, jac->measuretype);
375:     PetscCallExternal(HYPRE_BoomerAMGSetAggNumLevels, jac->hsolver, jac->agg_nl);
376:     PetscCallExternal(HYPRE_BoomerAMGSetPMaxElmts, jac->hsolver, jac->pmax);
377:     PetscCallExternal(HYPRE_BoomerAMGSetNumPaths, jac->hsolver, jac->agg_num_paths);
378:     PetscCallExternal(HYPRE_BoomerAMGSetCycleNumSweeps, jac->hsolver, jac->gridsweeps[0], 1);
379:     PetscCallExternal(HYPRE_BoomerAMGSetCycleNumSweeps, jac->hsolver, jac->gridsweeps[1], 2);
380:     PetscCallExternal(HYPRE_BoomerAMGSetCycleNumSweeps, jac->hsolver, jac->gridsweeps[2], 3);
381:     PetscCallExternal(HYPRE_BoomerAMGSetMaxCoarseSize, jac->hsolver, jac->maxc);
382:     PetscCallExternal(HYPRE_BoomerAMGSetMinCoarseSize, jac->hsolver, jac->minc);
383:     PetscCallExternal(HYPRE_BoomerAMGSetCoarsenType, jac->hsolver, jac->coarsentype);
384:     PetscCallExternal(HYPRE_BoomerAMGSetRelaxOrder, jac->hsolver, jac->relaxorder);
385:     PetscCallExternal(HYPRE_BoomerAMGSetInterpType, jac->hsolver, jac->interptype);
386:     PetscCallExternal(HYPRE_BoomerAMGSetRelaxType, jac->hsolver, jac->relaxtype[0]);
387:     /* GPU */
388: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
389:     PetscCall(PCMGGalerkinSetMatProductAlgorithm_HYPRE_BoomerAMG(pc, jac->spgemm_type));
390: #endif
391: #if PETSC_PKG_HYPRE_VERSION_GE(2, 18, 0)
392:     PetscCallExternal(HYPRE_BoomerAMGSetKeepTranspose, jac->hsolver, jac->keeptranspose ? 1 : 0);
393:     PetscCallExternal(HYPRE_BoomerAMGSetRAP2, jac->hsolver, jac->rap2);
394:     PetscCallExternal(HYPRE_BoomerAMGSetModuleRAP2, jac->hsolver, jac->mod_rap2);
395: #endif
396:     PetscCallExternal(HYPRE_BoomerAMGSetAggInterpType, jac->hsolver, jac->agg_interptype);

398:     /* AIR */
399: #if PETSC_PKG_HYPRE_VERSION_GE(2, 18, 0)
400:     PetscCallExternal(HYPRE_BoomerAMGSetRestriction, jac->hsolver, jac->Rtype);
401:     PetscCallExternal(HYPRE_BoomerAMGSetStrongThresholdR, jac->hsolver, jac->Rstrongthreshold);
402:     PetscCallExternal(HYPRE_BoomerAMGSetFilterThresholdR, jac->hsolver, jac->Rfilterthreshold);
403:     PetscCallExternal(HYPRE_BoomerAMGSetADropTol, jac->hsolver, jac->Adroptol);
404:     PetscCallExternal(HYPRE_BoomerAMGSetADropType, jac->hsolver, jac->Adroptype);
405: #endif

407:     PetscCall(MatGetBlockSize(pc->pmat, &bs));
408:     if (bs > 1) PetscCallExternal(HYPRE_BoomerAMGSetNumFunctions, jac->hsolver, bs);
409:     PetscCall(MatGetNearNullSpace(pc->mat, &mnull));
410:     if (mnull) {
411:       PetscCall(PCHYPREResetNearNullSpace_Private(pc));
412:       PetscCall(MatNullSpaceGetVecs(mnull, &has_const, &nvec, &vecs));
413:       PetscCall(PetscMalloc1(nvec + 1, &jac->hmnull));
414:       PetscCall(PetscMalloc1(nvec + 1, &jac->phmnull));
415:       for (i = 0; i < nvec; i++) {
416:         PetscCall(VecHYPRE_IJVectorCreate(vecs[i]->map, &jac->hmnull[i]));
417:         PetscCall(VecHYPRE_IJVectorCopy(vecs[i], jac->hmnull[i]));
418:         PetscCallExternal(HYPRE_IJVectorGetObject, jac->hmnull[i]->ij, (void **)&jac->phmnull[i]);
419:       }
420:       if (has_const) {
421:         PetscCall(MatCreateVecs(pc->pmat, &jac->hmnull_constant, NULL));
422:         PetscCall(VecSet(jac->hmnull_constant, 1));
423:         PetscCall(VecNormalize(jac->hmnull_constant, NULL));
424:         PetscCall(VecHYPRE_IJVectorCreate(jac->hmnull_constant->map, &jac->hmnull[nvec]));
425:         PetscCall(VecHYPRE_IJVectorCopy(jac->hmnull_constant, jac->hmnull[nvec]));
426:         PetscCallExternal(HYPRE_IJVectorGetObject, jac->hmnull[nvec]->ij, (void **)&jac->phmnull[nvec]);
427:         nvec++;
428:       }
429:       PetscCallExternal(HYPRE_BoomerAMGSetInterpVectors, jac->hsolver, nvec, jac->phmnull);
430:       jac->n_hmnull = nvec;
431:     }
432:   }

434:   /* special case for AMS */
435:   if (jac->setup == HYPRE_AMSSetup) {
436:     Mat_HYPRE         *hm;
437:     HYPRE_ParCSRMatrix parcsr;
438:     PetscCheck(jac->coords[0] || jac->constants[0] || jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1]), 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()");
439:     if (jac->dim) PetscCallExternal(HYPRE_AMSSetDimension, jac->hsolver, jac->dim);
440:     if (jac->constants[0]) {
441:       HYPRE_ParVector ozz, zoz, zzo = NULL;
442:       PetscCallExternal(HYPRE_IJVectorGetObject, jac->constants[0]->ij, (void **)(&ozz));
443:       PetscCallExternal(HYPRE_IJVectorGetObject, jac->constants[1]->ij, (void **)(&zoz));
444:       if (jac->constants[2]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->constants[2]->ij, (void **)(&zzo));
445:       PetscCallExternal(HYPRE_AMSSetEdgeConstantVectors, jac->hsolver, ozz, zoz, zzo);
446:     }
447:     if (jac->coords[0]) {
448:       HYPRE_ParVector coords[3];
449:       coords[0] = NULL;
450:       coords[1] = NULL;
451:       coords[2] = NULL;
452:       if (jac->coords[0]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[0]->ij, (void **)(&coords[0]));
453:       if (jac->coords[1]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[1]->ij, (void **)(&coords[1]));
454:       if (jac->coords[2]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[2]->ij, (void **)(&coords[2]));
455:       PetscCallExternal(HYPRE_AMSSetCoordinateVectors, jac->hsolver, coords[0], coords[1], coords[2]);
456:     }
457:     PetscCheck(jac->G, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "HYPRE AMS preconditioner needs the discrete gradient operator via PCHYPRESetDiscreteGradient");
458:     hm = (Mat_HYPRE *)jac->G->data;
459:     PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&parcsr));
460:     PetscCallExternal(HYPRE_AMSSetDiscreteGradient, jac->hsolver, parcsr);
461:     if (jac->alpha_Poisson) {
462:       hm = (Mat_HYPRE *)jac->alpha_Poisson->data;
463:       PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&parcsr));
464:       PetscCallExternal(HYPRE_AMSSetAlphaPoissonMatrix, jac->hsolver, parcsr);
465:     }
466:     if (jac->ams_beta_is_zero) {
467:       PetscCallExternal(HYPRE_AMSSetBetaPoissonMatrix, jac->hsolver, NULL);
468:     } else if (jac->beta_Poisson) {
469:       hm = (Mat_HYPRE *)jac->beta_Poisson->data;
470:       PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&parcsr));
471:       PetscCallExternal(HYPRE_AMSSetBetaPoissonMatrix, jac->hsolver, parcsr);
472:     } else if (jac->ams_beta_is_zero_part) {
473:       if (jac->interior) {
474:         HYPRE_ParVector interior = NULL;
475:         PetscCallExternal(HYPRE_IJVectorGetObject, jac->interior->ij, (void **)(&interior));
476:         PetscCallExternal(HYPRE_AMSSetInteriorNodes, jac->hsolver, interior);
477:       } else {
478:         jac->ams_beta_is_zero_part = PETSC_FALSE;
479:       }
480:     }
481:     if (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1])) {
482:       PetscInt           i;
483:       HYPRE_ParCSRMatrix nd_parcsrfull, nd_parcsr[3];
484:       if (jac->ND_PiFull) {
485:         hm = (Mat_HYPRE *)jac->ND_PiFull->data;
486:         PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&nd_parcsrfull));
487:       } else {
488:         nd_parcsrfull = NULL;
489:       }
490:       for (i = 0; i < 3; ++i) {
491:         if (jac->ND_Pi[i]) {
492:           hm = (Mat_HYPRE *)jac->ND_Pi[i]->data;
493:           PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&nd_parcsr[i]));
494:         } else {
495:           nd_parcsr[i] = NULL;
496:         }
497:       }
498:       PetscCallExternal(HYPRE_AMSSetInterpolations, jac->hsolver, nd_parcsrfull, nd_parcsr[0], nd_parcsr[1], nd_parcsr[2]);
499:     }
500:   }
501:   /* special case for ADS */
502:   if (jac->setup == HYPRE_ADSSetup) {
503:     Mat_HYPRE         *hm;
504:     HYPRE_ParCSRMatrix parcsr;
505:     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])))) {
506:       SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "HYPRE ADS preconditioner needs either the coordinate vectors via PCSetCoordinates() or the interpolation matrices via PCHYPRESetInterpolations");
507:     } 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");
508:     PetscCheck(jac->G, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "HYPRE ADS preconditioner needs the discrete gradient operator via PCHYPRESetDiscreteGradient");
509:     PetscCheck(jac->C, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "HYPRE ADS preconditioner needs the discrete curl operator via PCHYPRESetDiscreteGradient");
510:     if (jac->coords[0]) {
511:       HYPRE_ParVector coords[3];
512:       coords[0] = NULL;
513:       coords[1] = NULL;
514:       coords[2] = NULL;
515:       if (jac->coords[0]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[0]->ij, (void **)(&coords[0]));
516:       if (jac->coords[1]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[1]->ij, (void **)(&coords[1]));
517:       if (jac->coords[2]) PetscCallExternal(HYPRE_IJVectorGetObject, jac->coords[2]->ij, (void **)(&coords[2]));
518:       PetscCallExternal(HYPRE_ADSSetCoordinateVectors, jac->hsolver, coords[0], coords[1], coords[2]);
519:     }
520:     hm = (Mat_HYPRE *)jac->G->data;
521:     PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&parcsr));
522:     PetscCallExternal(HYPRE_ADSSetDiscreteGradient, jac->hsolver, parcsr);
523:     hm = (Mat_HYPRE *)jac->C->data;
524:     PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&parcsr));
525:     PetscCallExternal(HYPRE_ADSSetDiscreteCurl, jac->hsolver, parcsr);
526:     if ((jac->RT_PiFull || (jac->RT_Pi[0] && jac->RT_Pi[1])) && (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1]))) {
527:       PetscInt           i;
528:       HYPRE_ParCSRMatrix rt_parcsrfull, rt_parcsr[3];
529:       HYPRE_ParCSRMatrix nd_parcsrfull, nd_parcsr[3];
530:       if (jac->RT_PiFull) {
531:         hm = (Mat_HYPRE *)jac->RT_PiFull->data;
532:         PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&rt_parcsrfull));
533:       } else {
534:         rt_parcsrfull = NULL;
535:       }
536:       for (i = 0; i < 3; ++i) {
537:         if (jac->RT_Pi[i]) {
538:           hm = (Mat_HYPRE *)jac->RT_Pi[i]->data;
539:           PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&rt_parcsr[i]));
540:         } else {
541:           rt_parcsr[i] = NULL;
542:         }
543:       }
544:       if (jac->ND_PiFull) {
545:         hm = (Mat_HYPRE *)jac->ND_PiFull->data;
546:         PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&nd_parcsrfull));
547:       } else {
548:         nd_parcsrfull = NULL;
549:       }
550:       for (i = 0; i < 3; ++i) {
551:         if (jac->ND_Pi[i]) {
552:           hm = (Mat_HYPRE *)jac->ND_Pi[i]->data;
553:           PetscCallExternal(HYPRE_IJMatrixGetObject, hm->ij, (void **)(&nd_parcsr[i]));
554:         } else {
555:           nd_parcsr[i] = NULL;
556:         }
557:       }
558:       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]);
559:     }
560:   }
561:   PetscCallExternal(HYPRE_IJMatrixGetObject, hjac->ij, (void **)&hmat);
562:   PetscCallExternal(HYPRE_IJVectorGetObject, hjac->b->ij, (void **)&bv);
563:   PetscCallExternal(HYPRE_IJVectorGetObject, hjac->x->ij, (void **)&xv);
564:   PetscCall(PetscFPTrapPush(PETSC_FP_TRAP_OFF));
565:   PetscCallExternal(jac->setup, jac->hsolver, hmat, bv, xv);
566:   PetscCall(PetscFPTrapPop());
567:   PetscFunctionReturn(PETSC_SUCCESS);
568: }

570: static PetscErrorCode PCApply_HYPRE(PC pc, Vec b, Vec x)
571: {
572:   PC_HYPRE          *jac  = (PC_HYPRE *)pc->data;
573:   Mat_HYPRE         *hjac = (Mat_HYPRE *)jac->hpmat->data;
574:   HYPRE_ParCSRMatrix hmat;
575:   HYPRE_ParVector    jbv, jxv;

577:   PetscFunctionBegin;
578:   PetscCall(PetscCitationsRegister(hypreCitation, &cite));
579:   if (!jac->applyrichardson) PetscCall(VecSet(x, 0.0));
580:   PetscCall(VecHYPRE_IJVectorPushVecRead(hjac->b, b));
581:   if (jac->applyrichardson) PetscCall(VecHYPRE_IJVectorPushVec(hjac->x, x));
582:   else PetscCall(VecHYPRE_IJVectorPushVecWrite(hjac->x, x));
583:   PetscCallExternal(HYPRE_IJMatrixGetObject, hjac->ij, (void **)&hmat);
584:   PetscCallExternal(HYPRE_IJVectorGetObject, hjac->b->ij, (void **)&jbv);
585:   PetscCallExternal(HYPRE_IJVectorGetObject, hjac->x->ij, (void **)&jxv);
586:   PetscStackCallExternalVoid(
587:     "Hypre solve", do {
588:       HYPRE_Int hierr = (*jac->solve)(jac->hsolver, hmat, jbv, jxv);
589:       if (hierr) {
590:         PetscCheck(hierr == HYPRE_ERROR_CONV, PETSC_COMM_SELF, PETSC_ERR_LIB, "Error in HYPRE solver, error code %d", (int)hierr);
591:         HYPRE_ClearAllErrors();
592:       }
593:     } while (0));

595:   if (jac->setup == HYPRE_AMSSetup && jac->ams_beta_is_zero_part) PetscCallExternal(HYPRE_AMSProjectOutGradients, jac->hsolver, jxv);
596:   PetscCall(VecHYPRE_IJVectorPopVec(hjac->x));
597:   PetscCall(VecHYPRE_IJVectorPopVec(hjac->b));
598:   PetscFunctionReturn(PETSC_SUCCESS);
599: }

601: static PetscErrorCode PCMatApply_HYPRE_BoomerAMG(PC pc, Mat B, Mat X)
602: {
603:   PC_HYPRE           *jac  = (PC_HYPRE *)pc->data;
604:   Mat_HYPRE          *hjac = (Mat_HYPRE *)jac->hpmat->data;
605:   hypre_ParCSRMatrix *par_matrix;
606:   HYPRE_ParVector     hb, hx;
607:   const PetscScalar  *b;
608:   PetscScalar        *x;
609:   PetscInt            m, N, lda;
610:   hypre_Vector       *x_local;
611:   PetscMemType        type;

613:   PetscFunctionBegin;
614:   PetscCall(PetscCitationsRegister(hypreCitation, &cite));
615:   PetscCallExternal(HYPRE_IJMatrixGetObject, hjac->ij, (void **)&par_matrix);
616:   PetscCall(MatGetLocalSize(B, &m, NULL));
617:   PetscCall(MatGetSize(B, NULL, &N));
618:   PetscCallExternal(HYPRE_ParMultiVectorCreate, hypre_ParCSRMatrixComm(par_matrix), hypre_ParCSRMatrixGlobalNumRows(par_matrix), hypre_ParCSRMatrixRowStarts(par_matrix), N, &hb);
619:   PetscCallExternal(HYPRE_ParMultiVectorCreate, hypre_ParCSRMatrixComm(par_matrix), hypre_ParCSRMatrixGlobalNumRows(par_matrix), hypre_ParCSRMatrixRowStarts(par_matrix), N, &hx);
620:   PetscCall(MatZeroEntries(X));
621:   PetscCall(MatDenseGetArrayReadAndMemType(B, &b, &type));
622:   PetscCall(MatDenseGetLDA(B, &lda));
623:   PetscCheck(lda == m, PetscObjectComm((PetscObject)pc), PETSC_ERR_SUP, "Cannot use a LDA different than the number of local rows: % " PetscInt_FMT " != % " PetscInt_FMT, lda, m);
624:   PetscCall(MatDenseGetLDA(X, &lda));
625:   PetscCheck(lda == m, PetscObjectComm((PetscObject)pc), PETSC_ERR_SUP, "Cannot use a LDA different than the number of local rows: % " PetscInt_FMT " != % " PetscInt_FMT, lda, m);
626:   x_local = hypre_ParVectorLocalVector(hb);
627:   PetscCallExternal(hypre_SeqVectorSetDataOwner, x_local, 0);
628:   hypre_VectorData(x_local) = (HYPRE_Complex *)b;
629:   PetscCall(MatDenseGetArrayWriteAndMemType(X, &x, NULL));
630:   x_local = hypre_ParVectorLocalVector(hx);
631:   PetscCallExternal(hypre_SeqVectorSetDataOwner, x_local, 0);
632:   hypre_VectorData(x_local) = (HYPRE_Complex *)x;
633:   PetscCallExternal(hypre_ParVectorInitialize_v2, hb, type == PETSC_MEMTYPE_HOST ? HYPRE_MEMORY_HOST : HYPRE_MEMORY_DEVICE);
634:   PetscCallExternal(hypre_ParVectorInitialize_v2, hx, type == PETSC_MEMTYPE_HOST ? HYPRE_MEMORY_HOST : HYPRE_MEMORY_DEVICE);
635:   PetscStackCallExternalVoid(
636:     "Hypre solve", do {
637:       HYPRE_Int hierr = (*jac->solve)(jac->hsolver, par_matrix, hb, hx);
638:       if (hierr) {
639:         PetscCheck(hierr == HYPRE_ERROR_CONV, PETSC_COMM_SELF, PETSC_ERR_LIB, "Error in HYPRE solver, error code %d", (int)hierr);
640:         HYPRE_ClearAllErrors();
641:       }
642:     } while (0));
643:   PetscCallExternal(HYPRE_ParVectorDestroy, hb);
644:   PetscCallExternal(HYPRE_ParVectorDestroy, hx);
645:   PetscCall(MatDenseRestoreArrayReadAndMemType(B, &b));
646:   PetscCall(MatDenseRestoreArrayWriteAndMemType(X, &x));
647:   PetscFunctionReturn(PETSC_SUCCESS);
648: }

650: static PetscErrorCode PCReset_HYPRE(PC pc)
651: {
652:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;

654:   PetscFunctionBegin;
655:   PetscCall(MatDestroy(&jac->hpmat));
656:   PetscCall(MatDestroy(&jac->G));
657:   PetscCall(MatDestroy(&jac->C));
658:   PetscCall(MatDestroy(&jac->alpha_Poisson));
659:   PetscCall(MatDestroy(&jac->beta_Poisson));
660:   PetscCall(MatDestroy(&jac->RT_PiFull));
661:   PetscCall(MatDestroy(&jac->RT_Pi[0]));
662:   PetscCall(MatDestroy(&jac->RT_Pi[1]));
663:   PetscCall(MatDestroy(&jac->RT_Pi[2]));
664:   PetscCall(MatDestroy(&jac->ND_PiFull));
665:   PetscCall(MatDestroy(&jac->ND_Pi[0]));
666:   PetscCall(MatDestroy(&jac->ND_Pi[1]));
667:   PetscCall(MatDestroy(&jac->ND_Pi[2]));
668:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[0]));
669:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[1]));
670:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[2]));
671:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[0]));
672:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[1]));
673:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[2]));
674:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->interior));
675:   PetscCall(PCHYPREResetNearNullSpace_Private(pc));
676:   jac->ams_beta_is_zero      = PETSC_FALSE;
677:   jac->ams_beta_is_zero_part = PETSC_FALSE;
678:   jac->dim                   = 0;
679:   PetscFunctionReturn(PETSC_SUCCESS);
680: }

682: static PetscErrorCode PCDestroy_HYPRE(PC pc)
683: {
684:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;

686:   PetscFunctionBegin;
687:   PetscCall(PCReset_HYPRE(pc));
688:   if (jac->destroy) PetscCallExternal(jac->destroy, jac->hsolver);
689:   PetscCall(PetscFree(jac->hypre_type));
690:   if (jac->comm_hypre != MPI_COMM_NULL) PetscCall(PetscCommRestoreComm(PetscObjectComm((PetscObject)pc), &jac->comm_hypre));
691:   PetscCall(PetscFree(pc->data));

693:   PetscCall(PetscObjectChangeTypeName((PetscObject)pc, 0));
694:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetType_C", NULL));
695:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPREGetType_C", NULL));
696:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetDiscreteGradient_C", NULL));
697:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetDiscreteCurl_C", NULL));
698:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetInterpolations_C", NULL));
699:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetConstantEdgeVectors_C", NULL));
700:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetPoissonMatrix_C", NULL));
701:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetEdgeConstantVectors_C", NULL));
702:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPREAMSSetInteriorNodes_C", NULL));
703:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCGetInterpolations_C", NULL));
704:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCGetCoarseOperators_C", NULL));
705:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPREGetCFMarkers_C", NULL));
706:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCMGGalerkinSetMatProductAlgorithm_C", NULL));
707:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCMGGalerkinGetMatProductAlgorithm_C", NULL));
708:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSetCoordinates_C", NULL));
709:   PetscFunctionReturn(PETSC_SUCCESS);
710: }

712: static PetscErrorCode PCSetFromOptions_HYPRE_Pilut(PC pc, PetscOptionItems PetscOptionsObject)
713: {
714:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
715:   PetscBool flag;

717:   PetscFunctionBegin;
718:   PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE Pilut Options");
719:   PetscCall(PetscOptionsInt("-pc_hypre_pilut_maxiter", "Number of iterations", "None", jac->maxiter, &jac->maxiter, &flag));
720:   if (flag) PetscCallExternal(HYPRE_ParCSRPilutSetMaxIter, jac->hsolver, jac->maxiter);
721:   PetscCall(PetscOptionsReal("-pc_hypre_pilut_tol", "Drop tolerance", "None", jac->tol, &jac->tol, &flag));
722:   if (flag) PetscCallExternal(HYPRE_ParCSRPilutSetDropTolerance, jac->hsolver, jac->tol);
723:   PetscCall(PetscOptionsInt("-pc_hypre_pilut_factorrowsize", "FactorRowSize", "None", jac->factorrowsize, &jac->factorrowsize, &flag));
724:   if (flag) PetscCallExternal(HYPRE_ParCSRPilutSetFactorRowSize, jac->hsolver, jac->factorrowsize);
725:   PetscOptionsHeadEnd();
726:   PetscFunctionReturn(PETSC_SUCCESS);
727: }

729: static PetscErrorCode PCView_HYPRE_Pilut(PC pc, PetscViewer viewer)
730: {
731:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
732:   PetscBool isascii;

734:   PetscFunctionBegin;
735:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
736:   if (isascii) {
737:     PetscCall(PetscViewerASCIIPrintf(viewer, "  HYPRE Pilut preconditioning\n"));
738:     if (jac->maxiter != PETSC_DEFAULT) {
739:       PetscCall(PetscViewerASCIIPrintf(viewer, "    maximum number of iterations %" PetscInt_FMT "\n", jac->maxiter));
740:     } else {
741:       PetscCall(PetscViewerASCIIPrintf(viewer, "    default maximum number of iterations \n"));
742:     }
743:     if (jac->tol != PETSC_DEFAULT) {
744:       PetscCall(PetscViewerASCIIPrintf(viewer, "    drop tolerance %g\n", (double)jac->tol));
745:     } else {
746:       PetscCall(PetscViewerASCIIPrintf(viewer, "    default drop tolerance \n"));
747:     }
748:     if (jac->factorrowsize != PETSC_DEFAULT) {
749:       PetscCall(PetscViewerASCIIPrintf(viewer, "    factor row size %" PetscInt_FMT "\n", jac->factorrowsize));
750:     } else {
751:       PetscCall(PetscViewerASCIIPrintf(viewer, "    default factor row size \n"));
752:     }
753:   }
754:   PetscFunctionReturn(PETSC_SUCCESS);
755: }

757: static const char *HYPREILUType[] = {
758:   "Block-Jacobi-ILUk", "Block-Jacobi-ILUT", "", "", "", "", "", "", "", "", /* 0-9 */
759:   "GMRES-ILUk",        "GMRES-ILUT",        "", "", "", "", "", "", "", "", /* 10-19 */
760:   "NSH-ILUk",          "NSH-ILUT",          "", "", "", "", "", "", "", "", /* 20-29 */
761:   "RAS-ILUk",          "RAS-ILUT",          "", "", "", "", "", "", "", "", /* 30-39 */
762:   "ddPQ-GMRES-ILUk",   "ddPQ-GMRES-ILUT",   "", "", "", "", "", "", "", "", /* 40-49 */
763:   "GMRES-ILU0"                                                              /* 50 */
764: };

766: static const char *HYPREILUIterSetup[] = {"default", "async-in-place", "async-explicit", "sync-explicit", "semisync-explicit"};

768: static PetscErrorCode PCSetFromOptions_HYPRE_ILU(PC pc, PetscOptionItems PetscOptionsObject)
769: {
770:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
771:   PetscBool flg;
772:   PetscInt  indx;
773:   PetscReal tmpdbl;
774:   PetscBool tmp_truth;

776:   PetscFunctionBegin;
777:   PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE ILU Options");

779:   /* ILU: ILU Type */
780:   PetscCall(PetscOptionsEList("-pc_hypre_ilu_type", "Choose ILU Type", "None", HYPREILUType, PETSC_STATIC_ARRAY_LENGTH(HYPREILUType), HYPREILUType[0], &indx, &flg));
781:   if (flg) PetscCallExternal(HYPRE_ILUSetType, jac->hsolver, indx);

783:   /* ILU: ILU iterative setup type*/
784:   PetscCall(PetscOptionsEList("-pc_hypre_ilu_iterative_setup_type", "Set ILU iterative setup type", "None", HYPREILUIterSetup, PETSC_STATIC_ARRAY_LENGTH(HYPREILUIterSetup), HYPREILUIterSetup[0], &indx, &flg));
785:   if (flg) PetscCallExternal(HYPRE_ILUSetIterativeSetupType, jac->hsolver, indx);

787:   /* ILU: ILU iterative setup option*/
788:   PetscCall(PetscOptionsInt("-pc_hypre_ilu_iterative_setup_option", "Set ILU iterative setup option", "None", 0, &indx, &flg));
789:   if (flg) PetscCallExternal(HYPRE_ILUSetIterativeSetupOption, jac->hsolver, indx);

791:   /* ILU: ILU iterative setup maxiter */
792:   PetscCall(PetscOptionsInt("-pc_hypre_ilu_iterative_setup_maxiter", "Set ILU iterative setup maximum iteration count", "None", 0, &indx, &flg));
793:   if (flg) PetscCallExternal(HYPRE_ILUSetIterativeSetupMaxIter, jac->hsolver, indx);

795:   /* ILU: ILU iterative setup tolerance */
796:   PetscCall(PetscOptionsReal("-pc_hypre_ilu_iterative_setup_tolerance", "Set ILU iterative setup tolerance", "None", 0, &tmpdbl, &flg));
797:   if (flg) PetscCallExternal(HYPRE_ILUSetIterativeSetupTolerance, jac->hsolver, tmpdbl);

799:   /* ILU: ILU Print Level */
800:   PetscCall(PetscOptionsInt("-pc_hypre_ilu_print_level", "Set ILU print level", "None", 0, &indx, &flg));
801:   if (flg) PetscCallExternal(HYPRE_ILUSetPrintLevel, jac->hsolver, indx);

803:   /* ILU: Logging */
804:   PetscCall(PetscOptionsInt("-pc_hypre_ilu_logging", "Set ILU logging level", "None", 0, &indx, &flg));
805:   if (flg) PetscCallExternal(HYPRE_ILUSetLogging, jac->hsolver, indx);

807:   /* ILU: ILU Level */
808:   PetscCall(PetscOptionsInt("-pc_hypre_ilu_level", "Set ILU level", "None", 0, &indx, &flg));
809:   if (flg) PetscCallExternal(HYPRE_ILUSetLevelOfFill, jac->hsolver, indx);

811:   /* ILU: ILU Max NNZ per row */
812:   PetscCall(PetscOptionsInt("-pc_hypre_ilu_max_nnz_per_row", "Set maximum NNZ per row", "None", 0, &indx, &flg));
813:   if (flg) PetscCallExternal(HYPRE_ILUSetMaxNnzPerRow, jac->hsolver, indx);

815:   /* ILU: tolerance */
816:   PetscCall(PetscOptionsReal("-pc_hypre_ilu_tol", "Tolerance for ILU", "None", 0, &tmpdbl, &flg));
817:   if (flg) PetscCallExternal(HYPRE_ILUSetTol, jac->hsolver, tmpdbl);

819:   /* ILU: maximum iteration count */
820:   PetscCall(PetscOptionsInt("-pc_hypre_ilu_maxiter", "Set ILU max iterations", "None", 0, &indx, &flg));
821:   if (flg) PetscCallExternal(HYPRE_ILUSetMaxIter, jac->hsolver, indx);

823:   /* ILU: drop threshold */
824:   PetscCall(PetscOptionsReal("-pc_hypre_ilu_drop_threshold", "Drop threshold for ILU", "None", 0, &tmpdbl, &flg));
825:   if (flg) PetscCallExternal(HYPRE_ILUSetDropThreshold, jac->hsolver, tmpdbl);

827:   /* ILU: Triangular Solve */
828:   PetscCall(PetscOptionsBool("-pc_hypre_ilu_tri_solve", "Enable triangular solve", "None", PETSC_FALSE, &tmp_truth, &flg));
829:   if (flg) PetscCallExternal(HYPRE_ILUSetTriSolve, jac->hsolver, tmp_truth);

831:   /* ILU: Lower Jacobi iteration */
832:   PetscCall(PetscOptionsInt("-pc_hypre_ilu_lower_jacobi_iters", "Set lower Jacobi iteration count", "None", 0, &indx, &flg));
833:   if (flg) PetscCallExternal(HYPRE_ILUSetLowerJacobiIters, jac->hsolver, indx);

835:   /* ILU: Upper Jacobi iteration */
836:   PetscCall(PetscOptionsInt("-pc_hypre_ilu_upper_jacobi_iters", "Set upper Jacobi iteration count", "None", 0, &indx, &flg));
837:   if (flg) PetscCallExternal(HYPRE_ILUSetUpperJacobiIters, jac->hsolver, indx);

839:   /* ILU: local reordering */
840:   PetscCall(PetscOptionsBool("-pc_hypre_ilu_local_reordering", "Enable local reordering", "None", PETSC_FALSE, &tmp_truth, &flg));
841:   if (flg) PetscCallExternal(HYPRE_ILUSetLocalReordering, jac->hsolver, tmp_truth);

843:   PetscOptionsHeadEnd();
844:   PetscFunctionReturn(PETSC_SUCCESS);
845: }

847: static PetscErrorCode PCView_HYPRE_ILU(PC pc, PetscViewer viewer)
848: {
849:   PC_HYPRE         *jac      = (PC_HYPRE *)pc->data;
850:   hypre_ParILUData *ilu_data = (hypre_ParILUData *)jac->hsolver;
851:   PetscBool         isascii;
852:   PetscInt          indx;
853:   PetscReal         tmpdbl;
854:   PetscReal        *tmpdbl3;

856:   PetscFunctionBegin;
857:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
858:   if (isascii) {
859:     PetscCall(PetscViewerASCIIPrintf(viewer, "  HYPRE ILU preconditioning\n"));
860:     PetscStackCallExternalVoid("hypre_ParILUDataIluType", indx = hypre_ParILUDataIluType(ilu_data));
861:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU type              %s (%" PetscInt_FMT ")\n", HYPREILUType[indx], indx));
862:     PetscStackCallExternalVoid("hypre_ParILUDataLfil", indx = hypre_ParILUDataLfil(ilu_data));
863:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU level             %" PetscInt_FMT "\n", indx));
864:     PetscStackCallExternalVoid("hypre_ParILUDataMaxIter", indx = hypre_ParILUDataMaxIter(ilu_data));
865:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU max iterations    %" PetscInt_FMT "\n", indx));
866:     PetscStackCallExternalVoid("hypre_ParILUDataMaxRowNnz", indx = hypre_ParILUDataMaxRowNnz(ilu_data));
867:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU max NNZ per row   %" PetscInt_FMT "\n", indx));
868:     PetscStackCallExternalVoid("hypre_ParILUDataTriSolve", indx = hypre_ParILUDataTriSolve(ilu_data));
869:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU triangular solve  %" PetscInt_FMT "\n", indx));
870:     PetscStackCallExternalVoid("hypre_ParILUDataTol", tmpdbl = hypre_ParILUDataTol(ilu_data));
871:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU tolerance         %e\n", tmpdbl));
872:     PetscStackCallExternalVoid("hypre_ParILUDataDroptol", tmpdbl3 = hypre_ParILUDataDroptol(ilu_data));
873:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU drop tolerance    %e / %e / %e\n", tmpdbl3[0], tmpdbl3[1], tmpdbl3[2]));
874:     PetscStackCallExternalVoid("hypre_ParILUDataReorderingType", indx = hypre_ParILUDataReorderingType(ilu_data));
875:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU local reordering  %" PetscInt_FMT "\n", indx));
876:     PetscStackCallExternalVoid("hypre_ParILUDataLowerJacobiIters", indx = hypre_ParILUDataLowerJacobiIters(ilu_data));
877:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU lower Jacobi iterations  %" PetscInt_FMT "\n", indx));
878:     PetscStackCallExternalVoid("hypre_ParILUDataUpperJacobiIters", indx = hypre_ParILUDataUpperJacobiIters(ilu_data));
879:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU upper Jacobi iterations  %" PetscInt_FMT "\n", indx));
880:     PetscStackCallExternalVoid("hypre_ParILUDataPrintLevel", indx = hypre_ParILUDataPrintLevel(ilu_data));
881:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU print level      %" PetscInt_FMT "\n", indx));
882:     PetscStackCallExternalVoid("hypre_ParILUDataLogging", indx = hypre_ParILUDataLogging(ilu_data));
883:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU logging level    %" PetscInt_FMT "\n", indx));
884:     PetscStackCallExternalVoid("hypre_ParILUDataIterativeSetupType", indx = hypre_ParILUDataIterativeSetupType(ilu_data));
885:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU iterative setup type           %s (%" PetscInt_FMT ")\n", HYPREILUIterSetup[indx], indx));
886:     PetscStackCallExternalVoid("hypre_ParILUDataIterativeSetupOption", indx = hypre_ParILUDataIterativeSetupOption(ilu_data));
887:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU iterative setup option         %" PetscInt_FMT "\n", indx));
888:     PetscStackCallExternalVoid("hypre_ParILUDataIterativeSetupMaxIter", indx = hypre_ParILUDataIterativeSetupMaxIter(ilu_data));
889:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU iterative setup max iterations %" PetscInt_FMT "\n", indx));
890:     PetscStackCallExternalVoid("hypre_ParILUDataIterativeSetupTolerance", tmpdbl = hypre_ParILUDataIterativeSetupTolerance(ilu_data));
891:     PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU iterative setup tolerance      %e\n", tmpdbl));
892:   }
893:   PetscFunctionReturn(PETSC_SUCCESS);
894: }

896: static PetscErrorCode PCSetFromOptions_HYPRE_Euclid(PC pc, PetscOptionItems PetscOptionsObject)
897: {
898:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
899:   PetscBool flag, eu_bj = jac->eu_bj ? PETSC_TRUE : PETSC_FALSE;

901:   PetscFunctionBegin;
902:   PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE Euclid Options");
903:   PetscCall(PetscOptionsInt("-pc_hypre_euclid_level", "Factorization levels", "None", jac->eu_level, &jac->eu_level, &flag));
904:   if (flag) PetscCallExternal(HYPRE_EuclidSetLevel, jac->hsolver, jac->eu_level);

906:   PetscCall(PetscOptionsReal("-pc_hypre_euclid_droptolerance", "Drop tolerance for ILU(k) in Euclid", "None", jac->eu_droptolerance, &jac->eu_droptolerance, &flag));
907:   if (flag) {
908:     PetscMPIInt size;

910:     PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)pc), &size));
911:     PetscCheck(size == 1, PetscObjectComm((PetscObject)pc), PETSC_ERR_SUP, "hypre's Euclid does not support a parallel drop tolerance");
912:     PetscCallExternal(HYPRE_EuclidSetILUT, jac->hsolver, jac->eu_droptolerance);
913:   }

915:   PetscCall(PetscOptionsBool("-pc_hypre_euclid_bj", "Use Block Jacobi for ILU in Euclid", "None", eu_bj, &eu_bj, &flag));
916:   if (flag) {
917:     jac->eu_bj = eu_bj ? 1 : 0;
918:     PetscCallExternal(HYPRE_EuclidSetBJ, jac->hsolver, jac->eu_bj);
919:   }
920:   PetscOptionsHeadEnd();
921:   PetscFunctionReturn(PETSC_SUCCESS);
922: }

924: static PetscErrorCode PCView_HYPRE_Euclid(PC pc, PetscViewer viewer)
925: {
926:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
927:   PetscBool isascii;

929:   PetscFunctionBegin;
930:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
931:   if (isascii) {
932:     PetscCall(PetscViewerASCIIPrintf(viewer, "  HYPRE Euclid preconditioning\n"));
933:     if (jac->eu_level != PETSC_DEFAULT) {
934:       PetscCall(PetscViewerASCIIPrintf(viewer, "    factorization levels %" PetscInt_FMT "\n", jac->eu_level));
935:     } else {
936:       PetscCall(PetscViewerASCIIPrintf(viewer, "    default factorization levels \n"));
937:     }
938:     PetscCall(PetscViewerASCIIPrintf(viewer, "    drop tolerance %g\n", (double)jac->eu_droptolerance));
939:     PetscCall(PetscViewerASCIIPrintf(viewer, "    use Block-Jacobi? %" PetscInt_FMT "\n", jac->eu_bj));
940:   }
941:   PetscFunctionReturn(PETSC_SUCCESS);
942: }

944: static PetscErrorCode PCApplyTranspose_HYPRE_BoomerAMG(PC pc, Vec b, Vec x)
945: {
946:   PC_HYPRE          *jac  = (PC_HYPRE *)pc->data;
947:   Mat_HYPRE         *hjac = (Mat_HYPRE *)jac->hpmat->data;
948:   HYPRE_ParCSRMatrix hmat;
949:   HYPRE_ParVector    jbv, jxv;

951:   PetscFunctionBegin;
952:   PetscCall(PetscCitationsRegister(hypreCitation, &cite));
953:   PetscCall(VecSet(x, 0.0));
954:   PetscCall(VecHYPRE_IJVectorPushVecRead(hjac->b, b));
955:   PetscCall(VecHYPRE_IJVectorPushVecWrite(hjac->x, x));

957:   PetscCallExternal(HYPRE_IJMatrixGetObject, hjac->ij, (void **)&hmat);
958:   PetscCallExternal(HYPRE_IJVectorGetObject, hjac->b->ij, (void **)&jbv);
959:   PetscCallExternal(HYPRE_IJVectorGetObject, hjac->x->ij, (void **)&jxv);

961:   PetscStackCallExternalVoid(
962:     "Hypre Transpose solve", do {
963:       HYPRE_Int hierr = HYPRE_BoomerAMGSolveT(jac->hsolver, hmat, jbv, jxv);
964:       if (hierr) {
965:         /* error code of 1 in BoomerAMG merely means convergence not achieved */
966:         PetscCheck(hierr == 1, PETSC_COMM_SELF, PETSC_ERR_LIB, "Error in HYPRE solver, error code %d", (int)hierr);
967:         HYPRE_ClearAllErrors();
968:       }
969:     } while (0));

971:   PetscCall(VecHYPRE_IJVectorPopVec(hjac->x));
972:   PetscCall(VecHYPRE_IJVectorPopVec(hjac->b));
973:   PetscFunctionReturn(PETSC_SUCCESS);
974: }

976: static PetscErrorCode PCMGGalerkinGetMatProductAlgorithm_HYPRE_BoomerAMG(PC pc, const char *spgemm[])
977: {
978:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;

980:   PetscFunctionBegin;
982: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
983:   *spgemm = jac->spgemm_type;
984: #endif
985:   PetscFunctionReturn(PETSC_SUCCESS);
986: }

988: static const char *HYPREBoomerAMGCycleType[]   = {"", "V", "W"};
989: static const char *HYPREBoomerAMGCoarsenType[] = {"CLJP", "Ruge-Stueben", "", "modifiedRuge-Stueben", "", "", "Falgout", "", "PMIS", "", "HMIS"};
990: static const char *HYPREBoomerAMGMeasureType[] = {"local", "global"};
991: /* The following corresponds to HYPRE_BoomerAMGSetRelaxType which has many missing numbers in the enum */
992: static const char *HYPREBoomerAMGSmoothType[] = {"ILU", "Schwarz-smoothers", "Pilut", "ParaSails", "Euclid"};
993: 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"};
994: 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"};

996: static PetscErrorCode PCSetFromOptions_HYPRE_BoomerAMG(PC pc, PetscOptionItems PetscOptionsObject)
997: {
998:   PC_HYPRE   *jac = (PC_HYPRE *)pc->data;
999:   PetscInt    bs, n, indx, level;
1000:   PetscBool   flg, tmp_truth;
1001:   PetscReal   tmpdbl, twodbl[2];
1002:   const char *symtlist[] = {"nonsymmetric", "SPD", "nonsymmetric,SPD"};

1004:   PetscFunctionBegin;
1005:   PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE BoomerAMG Options");
1006:   PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_cycle_type", "Cycle type", "None", HYPREBoomerAMGCycleType + 1, 2, HYPREBoomerAMGCycleType[jac->cycletype], &indx, &flg));
1007:   if (flg) {
1008:     jac->cycletype = indx + 1;
1009:     PetscCallExternal(HYPRE_BoomerAMGSetCycleType, jac->hsolver, jac->cycletype);
1010:   }
1011:   PetscCall(PetscOptionsBoundedInt("-pc_hypre_boomeramg_max_levels", "Number of levels (of grids) allowed", "None", jac->maxlevels, &jac->maxlevels, &flg, 2));
1012:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetMaxLevels, jac->hsolver, jac->maxlevels);
1013:   PetscCall(PetscOptionsBoundedInt("-pc_hypre_boomeramg_max_iter", "Maximum iterations used PER hypre call", "None", jac->maxiter, &jac->maxiter, &flg, 1));
1014:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetMaxIter, jac->hsolver, jac->maxiter);
1015:   PetscCall(PetscOptionsBoundedReal("-pc_hypre_boomeramg_tol", "Convergence tolerance PER hypre call (0.0 = use a fixed number of iterations)", "None", jac->tol, &jac->tol, &flg, 0.0));
1016:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetTol, jac->hsolver, jac->tol);
1017:   bs = 1;
1018:   if (pc->pmat) PetscCall(MatGetBlockSize(pc->pmat, &bs));
1019:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_numfunctions", "Number of functions", "HYPRE_BoomerAMGSetNumFunctions", bs, &bs, &flg));
1020:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetNumFunctions, jac->hsolver, bs);

1022:   PetscCall(PetscOptionsBoundedReal("-pc_hypre_boomeramg_truncfactor", "Truncation factor for interpolation (0=no truncation)", "None", jac->truncfactor, &jac->truncfactor, &flg, 0.0));
1023:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetTruncFactor, jac->hsolver, jac->truncfactor);

1025:   PetscCall(PetscOptionsBoundedInt("-pc_hypre_boomeramg_P_max", "Max elements per row for interpolation operator (0=unlimited)", "None", jac->pmax, &jac->pmax, &flg, 0));
1026:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetPMaxElmts, jac->hsolver, jac->pmax);

1028:   PetscCall(PetscOptionsRangeInt("-pc_hypre_boomeramg_agg_nl", "Number of levels of aggressive coarsening", "None", jac->agg_nl, &jac->agg_nl, &flg, 0, jac->maxlevels));
1029:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetAggNumLevels, jac->hsolver, jac->agg_nl);

1031:   PetscCall(PetscOptionsBoundedInt("-pc_hypre_boomeramg_agg_num_paths", "Number of paths for aggressive coarsening", "None", jac->agg_num_paths, &jac->agg_num_paths, &flg, 1));
1032:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetNumPaths, jac->hsolver, jac->agg_num_paths);

1034:   PetscCall(PetscOptionsBoundedReal("-pc_hypre_boomeramg_strong_threshold", "Threshold for being strongly connected", "None", jac->strongthreshold, &jac->strongthreshold, &flg, 0.0));
1035:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetStrongThreshold, jac->hsolver, jac->strongthreshold);
1036:   PetscCall(PetscOptionsRangeReal("-pc_hypre_boomeramg_max_row_sum", "Maximum row sum", "None", jac->maxrowsum, &jac->maxrowsum, &flg, 0.0, 1.0));
1037:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetMaxRowSum, jac->hsolver, jac->maxrowsum);

1039:   /* Grid sweeps */
1040:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_all", "Number of sweeps for the up and down grid levels", "None", jac->gridsweeps[0], &indx, &flg));
1041:   if (flg) {
1042:     /* modify the jac structure so we can view the updated options with PC_View */
1043:     jac->gridsweeps[0] = indx;
1044:     jac->gridsweeps[1] = indx;
1045:     /*defaults coarse to 1 */
1046:     jac->gridsweeps[2] = 1;
1047:   }
1048:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_nodal_coarsen", "Use a nodal based coarsening 1-6", "HYPRE_BoomerAMGSetNodal", jac->nodal_coarsening, &jac->nodal_coarsening, &flg));
1049:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetNodal, jac->hsolver, jac->nodal_coarsening);
1050:   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));
1051:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetNodalDiag, jac->hsolver, jac->nodal_coarsening_diag);
1052:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_vec_interp_variant", "Variant of algorithm 1-3", "HYPRE_BoomerAMGSetInterpVecVariant", jac->vec_interp_variant, &jac->vec_interp_variant, &flg));
1053:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetInterpVecVariant, jac->hsolver, jac->vec_interp_variant);
1054:   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));
1055:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetInterpVecQMax, jac->hsolver, jac->vec_interp_qmax);
1056:   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));
1057:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetSmoothInterpVectors, jac->hsolver, jac->vec_interp_smooth);
1058:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_interp_refine", "Preprocess the interpolation matrix through iterative weight refinement", "HYPRE_BoomerAMGSetInterpRefine", jac->interp_refine, &jac->interp_refine, &flg));
1059:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetInterpRefine, jac->hsolver, jac->interp_refine);
1060:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_down", "Number of sweeps for the down cycles", "None", jac->gridsweeps[0], &indx, &flg));
1061:   if (flg) {
1062:     PetscCallExternal(HYPRE_BoomerAMGSetCycleNumSweeps, jac->hsolver, indx, 1);
1063:     jac->gridsweeps[0] = indx;
1064:   }
1065:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_up", "Number of sweeps for the up cycles", "None", jac->gridsweeps[1], &indx, &flg));
1066:   if (flg) {
1067:     PetscCallExternal(HYPRE_BoomerAMGSetCycleNumSweeps, jac->hsolver, indx, 2);
1068:     jac->gridsweeps[1] = indx;
1069:   }
1070:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_coarse", "Number of sweeps for the coarse level", "None", jac->gridsweeps[2], &indx, &flg));
1071:   if (flg) {
1072:     PetscCallExternal(HYPRE_BoomerAMGSetCycleNumSweeps, jac->hsolver, indx, 3);
1073:     jac->gridsweeps[2] = indx;
1074:   }

1076:   /* Smooth type */
1077:   PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_smooth_type", "Enable more complex smoothers", "None", HYPREBoomerAMGSmoothType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGSmoothType), HYPREBoomerAMGSmoothType[0], &indx, &flg));
1078:   if (flg) {
1079:     jac->smoothtype = indx;
1080:     PetscCallExternal(HYPRE_BoomerAMGSetSmoothType, jac->hsolver, indx + 5);
1081:     jac->smoothnumlevels = 25;
1082:     PetscCallExternal(HYPRE_BoomerAMGSetSmoothNumLevels, jac->hsolver, 25);
1083:   }

1085:   /* Number of smoothing levels */
1086:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_smooth_num_levels", "Number of levels on which more complex smoothers are used", "None", 25, &indx, &flg));
1087:   if (flg && (jac->smoothtype != -1)) {
1088:     jac->smoothnumlevels = indx;
1089:     PetscCallExternal(HYPRE_BoomerAMGSetSmoothNumLevels, jac->hsolver, indx);
1090:   }

1092:   /* Smooth num sweeps */
1093:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_smooth_num_sweeps", "Set number of smoother sweeps", "None", 1, &indx, &flg));
1094:   if (flg && indx > 0) {
1095:     jac->smoothsweeps = indx;
1096:     PetscCallExternal(HYPRE_BoomerAMGSetSmoothNumSweeps, jac->hsolver, indx);
1097:   }

1099:   /* ILU: ILU Type */
1100:   PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_ilu_type", "Choose ILU Type", "None", HYPREILUType, PETSC_STATIC_ARRAY_LENGTH(HYPREILUType), HYPREILUType[0], &indx, &flg));
1101:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILUType, jac->hsolver, indx);

1103:   /* ILU: ILU iterative setup type*/
1104:   PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_ilu_iterative_setup_type", "Set ILU iterative setup type", "None", HYPREILUIterSetup, PETSC_STATIC_ARRAY_LENGTH(HYPREILUIterSetup), HYPREILUIterSetup[0], &indx, &flg));
1105:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILUIterSetupType, jac->hsolver, indx);

1107:   /* ILU: ILU iterative setup option*/
1108:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_ilu_iterative_setup_option", "Set ILU iterative setup option", "None", 0, &indx, &flg));
1109:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILUIterSetupOption, jac->hsolver, indx);

1111:   /* ILU: ILU iterative setup maxiter */
1112:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_ilu_iterative_setup_maxiter", "Set ILU iterative setup maximum iteration count", "None", 0, &indx, &flg));
1113:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILUIterSetupMaxIter, jac->hsolver, indx);

1115:   /* ILU: ILU iterative setup tolerance */
1116:   PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_ilu_iterative_setup_tolerance", "Set ILU iterative setup tolerance", "None", 0, &tmpdbl, &flg));
1117:   if (flg) PetscCallExternal(hypre_BoomerAMGSetILUIterSetupTolerance, jac->hsolver, tmpdbl);

1119:   /* ILU: ILU Print Level */
1120:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_ilu_print_level", "Set ILU print level", "None", 0, &indx, &flg));
1121:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetPrintLevel, jac->hsolver, indx);

1123:   /* ILU: Logging */
1124:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_ilu_logging", "Set ILU logging level", "None", 0, &indx, &flg));
1125:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetLogging, jac->hsolver, indx);

1127:   /* ILU: ILU Level */
1128:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_ilu_level", "Set ILU level", "None", 0, &indx, &flg));
1129:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILULevel, jac->hsolver, indx);

1131:   /* ILU: ILU Max NNZ per row */
1132:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_ilu_max_nnz_per_row", "Set maximum NNZ per row", "None", 0, &indx, &flg));
1133:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILUMaxRowNnz, jac->hsolver, indx);

1135:   /* ILU: maximum iteration count */
1136:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_ilu_maxiter", "Set ILU max iterations", "None", 0, &indx, &flg));
1137:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILUMaxIter, jac->hsolver, indx);

1139:   /* ILU: drop threshold */
1140:   PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_ilu_drop_tol", "Drop tolerance for ILU", "None", 0, &tmpdbl, &flg));
1141:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILUDroptol, jac->hsolver, tmpdbl);

1143:   /* ILU: Triangular Solve */
1144:   PetscCall(PetscOptionsBool("-pc_hypre_boomeramg_ilu_tri_solve", "Enable triangular solve", "None", PETSC_FALSE, &tmp_truth, &flg));
1145:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILUTriSolve, jac->hsolver, tmp_truth);

1147:   /* ILU: Lower Jacobi iteration */
1148:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_ilu_lower_jacobi_iters", "Set lower Jacobi iteration count", "None", 0, &indx, &flg));
1149:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILULowerJacobiIters, jac->hsolver, indx);

1151:   /* ILU: Upper Jacobi iteration */
1152:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_ilu_upper_jacobi_iters", "Set upper Jacobi iteration count", "None", 0, &indx, &flg));
1153:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILUUpperJacobiIters, jac->hsolver, indx);

1155:   /* ILU: local reordering */
1156:   PetscCall(PetscOptionsBool("-pc_hypre_boomeramg_ilu_local_reordering", "Enable local reordering", "None", PETSC_FALSE, &tmp_truth, &flg));
1157:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetILULocalReordering, jac->hsolver, tmp_truth);

1159:   /* Number of levels for ILU(k) for Euclid */
1160:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_eu_level", "Number of levels for ILU(k) in Euclid smoother", "None", 0, &indx, &flg));
1161:   if (flg && (jac->smoothtype == 4)) {
1162:     jac->eu_level = indx;
1163:     PetscCallExternal(HYPRE_BoomerAMGSetEuLevel, jac->hsolver, indx);
1164:   }

1166:   /* Filter for ILU(k) for Euclid */
1167:   PetscReal droptolerance;
1168:   PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_eu_droptolerance", "Drop tolerance for ILU(k) in Euclid smoother", "None", 0, &droptolerance, &flg));
1169:   if (flg && (jac->smoothtype == 4)) {
1170:     jac->eu_droptolerance = droptolerance;
1171:     PetscCallExternal(HYPRE_BoomerAMGSetEuLevel, jac->hsolver, droptolerance);
1172:   }

1174:   /* Use Block Jacobi ILUT for Euclid */
1175:   PetscCall(PetscOptionsBool("-pc_hypre_boomeramg_eu_bj", "Use Block Jacobi for ILU in Euclid smoother?", "None", PETSC_FALSE, &tmp_truth, &flg));
1176:   if (flg && (jac->smoothtype == 4)) {
1177:     jac->eu_bj = tmp_truth;
1178:     PetscCallExternal(HYPRE_BoomerAMGSetEuBJ, jac->hsolver, jac->eu_bj);
1179:   }

1181:   /* Relax type */
1182:   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));
1183:   if (flg) {
1184:     jac->relaxtype[0] = jac->relaxtype[1] = indx;
1185:     PetscCallExternal(HYPRE_BoomerAMGSetRelaxType, jac->hsolver, indx);
1186:     /* by default, coarse type set to 9 */
1187:     jac->relaxtype[2] = 9;
1188:     PetscCallExternal(HYPRE_BoomerAMGSetCycleRelaxType, jac->hsolver, 9, 3);
1189:   }
1190:   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));
1191:   if (flg) {
1192:     jac->relaxtype[0] = indx;
1193:     PetscCallExternal(HYPRE_BoomerAMGSetCycleRelaxType, jac->hsolver, indx, 1);
1194:   }
1195:   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));
1196:   if (flg) {
1197:     jac->relaxtype[1] = indx;
1198:     PetscCallExternal(HYPRE_BoomerAMGSetCycleRelaxType, jac->hsolver, indx, 2);
1199:   }
1200:   PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_relax_type_coarse", "Relax type on coarse grid", "None", HYPREBoomerAMGRelaxType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGRelaxType), HYPREBoomerAMGRelaxType[9], &indx, &flg));
1201:   if (flg) {
1202:     jac->relaxtype[2] = indx;
1203:     PetscCallExternal(HYPRE_BoomerAMGSetCycleRelaxType, jac->hsolver, indx, 3);
1204:   }

1206:   /* Relaxation Weight */
1207:   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));
1208:   if (flg) {
1209:     PetscCallExternal(HYPRE_BoomerAMGSetRelaxWt, jac->hsolver, tmpdbl);
1210:     jac->relaxweight = tmpdbl;
1211:   }

1213:   n         = 2;
1214:   twodbl[0] = twodbl[1] = 1.0;
1215:   PetscCall(PetscOptionsRealArray("-pc_hypre_boomeramg_relax_weight_level", "Set the relaxation weight for a particular level (weight,level)", "None", twodbl, &n, &flg));
1216:   if (flg) {
1217:     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);
1218:     indx = (int)PetscAbsReal(twodbl[1]);
1219:     PetscCallExternal(HYPRE_BoomerAMGSetLevelRelaxWt, jac->hsolver, twodbl[0], indx);
1220:   }

1222:   /* Outer relaxation Weight */
1223:   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));
1224:   if (flg) {
1225:     PetscCallExternal(HYPRE_BoomerAMGSetOuterWt, jac->hsolver, tmpdbl);
1226:     jac->outerrelaxweight = tmpdbl;
1227:   }

1229:   n         = 2;
1230:   twodbl[0] = twodbl[1] = 1.0;
1231:   PetscCall(PetscOptionsRealArray("-pc_hypre_boomeramg_outer_relax_weight_level", "Set the outer relaxation weight for a particular level (weight,level)", "None", twodbl, &n, &flg));
1232:   if (flg) {
1233:     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);
1234:     indx = (int)PetscAbsReal(twodbl[1]);
1235:     PetscCallExternal(HYPRE_BoomerAMGSetLevelOuterWt, jac->hsolver, twodbl[0], indx);
1236:   }

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

1241:   if (flg && tmp_truth) {
1242:     jac->relaxorder = 0;
1243:     PetscCallExternal(HYPRE_BoomerAMGSetRelaxOrder, jac->hsolver, jac->relaxorder);
1244:   }
1245:   PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_measure_type", "Measure type", "None", HYPREBoomerAMGMeasureType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGMeasureType), HYPREBoomerAMGMeasureType[0], &indx, &flg));
1246:   if (flg) {
1247:     jac->measuretype = indx;
1248:     PetscCallExternal(HYPRE_BoomerAMGSetMeasureType, jac->hsolver, jac->measuretype);
1249:   }
1250:   /* update list length 3/07 */
1251:   PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_coarsen_type", "Coarsen type", "None", HYPREBoomerAMGCoarsenType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGCoarsenType), HYPREBoomerAMGCoarsenType[6], &indx, &flg));
1252:   if (flg) {
1253:     jac->coarsentype = indx;
1254:     PetscCallExternal(HYPRE_BoomerAMGSetCoarsenType, jac->hsolver, jac->coarsentype);
1255:   }

1257:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_max_coarse_size", "Maximum size of coarsest grid", "None", jac->maxc, &jac->maxc, &flg));
1258:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetMaxCoarseSize, jac->hsolver, jac->maxc);
1259:   PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_min_coarse_size", "Minimum size of coarsest grid", "None", jac->minc, &jac->minc, &flg));
1260:   if (flg) PetscCallExternal(HYPRE_BoomerAMGSetMinCoarseSize, jac->hsolver, jac->minc);
1261: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
1262:   // global parameter but is closely associated with BoomerAMG
1263:   PetscCall(PetscOptionsEList("-pc_mg_galerkin_mat_product_algorithm", "Type of SpGEMM to use in hypre (only for now)", "PCMGGalerkinSetMatProductAlgorithm", HYPRESpgemmTypes, PETSC_STATIC_ARRAY_LENGTH(HYPRESpgemmTypes), HYPRESpgemmTypes[0], &indx, &flg));
1264:   if (flg) PetscCall(PCMGGalerkinSetMatProductAlgorithm_HYPRE_BoomerAMG(pc, HYPRESpgemmTypes[indx]));
1265: #endif
1266:   /* AIR */
1267: #if PETSC_PKG_HYPRE_VERSION_GE(2, 18, 0)
1268:   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));
1269:   PetscCallExternal(HYPRE_BoomerAMGSetRestriction, jac->hsolver, jac->Rtype);
1270:   if (jac->Rtype) {
1271:     HYPRE_Int **grid_relax_points = hypre_TAlloc(HYPRE_Int *, 4, HYPRE_MEMORY_HOST);
1272:     char       *prerelax[256];
1273:     char       *postrelax[256];
1274:     char        stringF[2] = "F", stringC[2] = "C", stringA[2] = "A";
1275:     PetscInt    ns_down = 256, ns_up = 256;
1276:     PetscBool   matchF, matchC, matchA;

1278:     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 */

1280:     PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_strongthresholdR", "Threshold for R", "None", jac->Rstrongthreshold, &jac->Rstrongthreshold, NULL));
1281:     PetscCallExternal(HYPRE_BoomerAMGSetStrongThresholdR, jac->hsolver, jac->Rstrongthreshold);

1283:     PetscCall(PetscOptionsReal("-pc_hypre_boomeramg_filterthresholdR", "Filter threshold for R", "None", jac->Rfilterthreshold, &jac->Rfilterthreshold, NULL));
1284:     PetscCallExternal(HYPRE_BoomerAMGSetFilterThresholdR, jac->hsolver, jac->Rfilterthreshold);

1286:     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));
1287:     PetscCallExternal(HYPRE_BoomerAMGSetADropTol, jac->hsolver, jac->Adroptol);

1289:     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));
1290:     PetscCallExternal(HYPRE_BoomerAMGSetADropType, jac->hsolver, jac->Adroptype);
1291:     PetscCall(PetscOptionsStringArray("-pc_hypre_boomeramg_prerelax", "Defines prerelax scheme", "None", prerelax, &ns_down, NULL));
1292:     PetscCall(PetscOptionsStringArray("-pc_hypre_boomeramg_postrelax", "Defines postrelax scheme", "None", postrelax, &ns_up, NULL));
1293:     PetscCheck(ns_down == jac->gridsweeps[0], PetscObjectComm((PetscObject)jac), PETSC_ERR_ARG_SIZ, "The number of arguments passed to -pc_hypre_boomeramg_prerelax must match the number passed to -pc_hypre_bomeramg_grid_sweeps_down");
1294:     PetscCheck(ns_up == jac->gridsweeps[1], PetscObjectComm((PetscObject)jac), PETSC_ERR_ARG_SIZ, "The number of arguments passed to -pc_hypre_boomeramg_postrelax must match the number passed to -pc_hypre_bomeramg_grid_sweeps_up");

1296:     grid_relax_points[0]    = NULL;
1297:     grid_relax_points[1]    = hypre_TAlloc(HYPRE_Int, ns_down, HYPRE_MEMORY_HOST);
1298:     grid_relax_points[2]    = hypre_TAlloc(HYPRE_Int, ns_up, HYPRE_MEMORY_HOST);
1299:     grid_relax_points[3]    = hypre_TAlloc(HYPRE_Int, jac->gridsweeps[2], HYPRE_MEMORY_HOST);
1300:     grid_relax_points[3][0] = 0;

1302:     // set down relax scheme
1303:     for (PetscInt i = 0; i < ns_down; i++) {
1304:       PetscCall(PetscStrcasecmp(prerelax[i], stringF, &matchF));
1305:       PetscCall(PetscStrcasecmp(prerelax[i], stringC, &matchC));
1306:       PetscCall(PetscStrcasecmp(prerelax[i], stringA, &matchA));
1307:       PetscCheck(matchF || matchC || matchA, PetscObjectComm((PetscObject)jac), PETSC_ERR_ARG_WRONG, "Valid argument options for -pc_hypre_boomeramg_prerelax are C, F, and A");
1308:       if (matchF) grid_relax_points[1][i] = -1;
1309:       else if (matchC) grid_relax_points[1][i] = 1;
1310:       else if (matchA) grid_relax_points[1][i] = 0;
1311:     }

1313:     // set up relax scheme
1314:     for (PetscInt i = 0; i < ns_up; i++) {
1315:       PetscCall(PetscStrcasecmp(postrelax[i], stringF, &matchF));
1316:       PetscCall(PetscStrcasecmp(postrelax[i], stringC, &matchC));
1317:       PetscCall(PetscStrcasecmp(postrelax[i], stringA, &matchA));
1318:       PetscCheck(matchF || matchC || matchA, PetscObjectComm((PetscObject)jac), PETSC_ERR_ARG_WRONG, "Valid argument options for -pc_hypre_boomeramg_postrelax are C, F, and A");
1319:       if (matchF) grid_relax_points[2][i] = -1;
1320:       else if (matchC) grid_relax_points[2][i] = 1;
1321:       else if (matchA) grid_relax_points[2][i] = 0;
1322:     }

1324:     // set coarse relax scheme
1325:     for (PetscInt i = 0; i < jac->gridsweeps[2]; i++) grid_relax_points[3][i] = 0;

1327:     // Pass relax schemes to hypre
1328:     PetscCallExternal(HYPRE_BoomerAMGSetGridRelaxPoints, jac->hsolver, grid_relax_points);

1330:     // cleanup memory
1331:     for (PetscInt i = 0; i < ns_down; i++) PetscCall(PetscFree(prerelax[i]));
1332:     for (PetscInt i = 0; i < ns_up; i++) PetscCall(PetscFree(postrelax[i]));
1333:   }
1334: #endif

1336: #if PETSC_PKG_HYPRE_VERSION_LE(9, 9, 9)
1337:   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);
1338: #endif

1340:   /* new 3/07 */
1341:   PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_interp_type", "Interpolation type", "None", HYPREBoomerAMGInterpType, PETSC_STATIC_ARRAY_LENGTH(HYPREBoomerAMGInterpType), HYPREBoomerAMGInterpType[0], &indx, &flg));
1342:   if (flg || jac->Rtype) {
1343:     if (flg) jac->interptype = indx;
1344:     PetscCallExternal(HYPRE_BoomerAMGSetInterpType, jac->hsolver, jac->interptype);
1345:   }

1347:   PetscCall(PetscOptionsName("-pc_hypre_boomeramg_print_statistics", "Print statistics", "None", &flg));
1348:   if (flg) {
1349:     level = 3;
1350:     PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_print_statistics", "Print statistics", "None", level, &level, NULL));

1352:     jac->printstatistics = PETSC_TRUE;
1353:     PetscCallExternal(HYPRE_BoomerAMGSetPrintLevel, jac->hsolver, level);
1354:   }

1356:   PetscCall(PetscOptionsName("-pc_hypre_boomeramg_print_debug", "Print debug information", "None", &flg));
1357:   if (flg) {
1358:     level = 3;
1359:     PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_print_debug", "Print debug information", "None", level, &level, NULL));

1361:     jac->printstatistics = PETSC_TRUE;
1362:     PetscCallExternal(HYPRE_BoomerAMGSetDebugFlag, jac->hsolver, level);
1363:   }

1365:   PetscCall(PetscOptionsBool("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None", PETSC_FALSE, &tmp_truth, &flg));
1366:   if (flg && tmp_truth) {
1367:     PetscInt tmp_int;
1368:     PetscCall(PetscOptionsInt("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None", jac->nodal_relax_levels, &tmp_int, &flg));
1369:     if (flg) jac->nodal_relax_levels = tmp_int;
1370:     PetscCallExternal(HYPRE_BoomerAMGSetSmoothType, jac->hsolver, 6);
1371:     PetscCallExternal(HYPRE_BoomerAMGSetDomainType, jac->hsolver, 1);
1372:     PetscCallExternal(HYPRE_BoomerAMGSetOverlap, jac->hsolver, 0);
1373:     PetscCallExternal(HYPRE_BoomerAMGSetSmoothNumLevels, jac->hsolver, jac->nodal_relax_levels);
1374:   }

1376:   PetscCall(PetscOptionsBool("-pc_hypre_boomeramg_keeptranspose", "Avoid transpose matvecs in preconditioner application", "None", jac->keeptranspose, &jac->keeptranspose, NULL));
1377:   PetscCallExternal(HYPRE_BoomerAMGSetKeepTranspose, jac->hsolver, jac->keeptranspose ? 1 : 0);

1379:   /* options for ParaSails solvers */
1380:   PetscCall(PetscOptionsEList("-pc_hypre_boomeramg_parasails_sym", "Symmetry of matrix and preconditioner", "None", symtlist, PETSC_STATIC_ARRAY_LENGTH(symtlist), symtlist[0], &indx, &flg));
1381:   if (flg) {
1382:     jac->symt = indx;
1383:     PetscCallExternal(HYPRE_BoomerAMGSetSym, jac->hsolver, jac->symt);
1384:   }

1386:   PetscOptionsHeadEnd();
1387:   PetscFunctionReturn(PETSC_SUCCESS);
1388: }

1390: 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)
1391: {
1392:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1393:   HYPRE_Int oits;

1395:   PetscFunctionBegin;
1396:   PetscCall(PetscCitationsRegister(hypreCitation, &cite));
1397:   PetscCallExternal(HYPRE_BoomerAMGSetMaxIter, jac->hsolver, its * jac->maxiter);
1398:   PetscCallExternal(HYPRE_BoomerAMGSetTol, jac->hsolver, rtol);
1399:   jac->applyrichardson = PETSC_TRUE;
1400:   PetscCall(PCApply_HYPRE(pc, b, y));
1401:   jac->applyrichardson = PETSC_FALSE;
1402:   PetscCallExternal(HYPRE_BoomerAMGGetNumIterations, jac->hsolver, &oits);
1403:   *outits = oits;
1404:   if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
1405:   else *reason = PCRICHARDSON_CONVERGED_RTOL;
1406:   PetscCallExternal(HYPRE_BoomerAMGSetTol, jac->hsolver, jac->tol);
1407:   PetscCallExternal(HYPRE_BoomerAMGSetMaxIter, jac->hsolver, jac->maxiter);
1408:   PetscFunctionReturn(PETSC_SUCCESS);
1409: }

1411: static PetscErrorCode PCView_HYPRE_BoomerAMG(PC pc, PetscViewer viewer)
1412: {
1413:   PC_HYPRE         *jac      = (PC_HYPRE *)pc->data;
1414:   hypre_ParAMGData *amg_data = (hypre_ParAMGData *)jac->hsolver;
1415:   PetscBool         isascii;
1416:   PetscInt          indx;
1417:   PetscReal         val;

1419:   PetscFunctionBegin;
1420:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
1421:   if (isascii) {
1422:     PetscCall(PetscViewerASCIIPrintf(viewer, "  HYPRE BoomerAMG preconditioning\n"));
1423:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Cycle type %s\n", HYPREBoomerAMGCycleType[jac->cycletype]));
1424:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Maximum number of levels %" PetscInt_FMT "\n", jac->maxlevels));
1425:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Maximum number of iterations PER hypre call %" PetscInt_FMT "\n", jac->maxiter));
1426:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Convergence tolerance PER hypre call %g\n", (double)jac->tol));
1427:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Threshold for strong coupling %g\n", (double)jac->strongthreshold));
1428:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Interpolation truncation factor %g\n", (double)jac->truncfactor));
1429:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Interpolation: max elements per row %" PetscInt_FMT "\n", jac->pmax));
1430:     if (jac->interp_refine) PetscCall(PetscViewerASCIIPrintf(viewer, "    Interpolation: number of steps of weighted refinement %" PetscInt_FMT "\n", jac->interp_refine));
1431:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Number of levels of aggressive coarsening %" PetscInt_FMT "\n", jac->agg_nl));
1432:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Number of paths for aggressive coarsening %" PetscInt_FMT "\n", jac->agg_num_paths));

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

1436:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Sweeps down         %" PetscInt_FMT "\n", jac->gridsweeps[0]));
1437:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Sweeps up           %" PetscInt_FMT "\n", jac->gridsweeps[1]));
1438:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Sweeps on coarse    %" PetscInt_FMT "\n", jac->gridsweeps[2]));

1440:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Relax down          %s\n", HYPREBoomerAMGRelaxType[jac->relaxtype[0]]));
1441:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Relax up            %s\n", HYPREBoomerAMGRelaxType[jac->relaxtype[1]]));
1442:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Relax on coarse     %s\n", HYPREBoomerAMGRelaxType[jac->relaxtype[2]]));

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

1447:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Maximum size of coarsest grid %" PetscInt_FMT "\n", jac->maxc));
1448:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Minimum size of coarsest grid %" PetscInt_FMT "\n", jac->minc));

1450:     if (jac->relaxorder) {
1451:       PetscCall(PetscViewerASCIIPrintf(viewer, "    Using CF-relaxation\n"));
1452:     } else {
1453:       PetscCall(PetscViewerASCIIPrintf(viewer, "    Not using CF-relaxation\n"));
1454:     }
1455:     if (jac->smoothtype != -1) {
1456:       PetscCall(PetscViewerASCIIPrintf(viewer, "    Smooth type          %s\n", HYPREBoomerAMGSmoothType[jac->smoothtype]));
1457:       PetscCall(PetscViewerASCIIPrintf(viewer, "    Smooth num levels    %" PetscInt_FMT "\n", jac->smoothnumlevels));
1458:       PetscCall(PetscViewerASCIIPrintf(viewer, "    Smooth num sweeps    %" PetscInt_FMT "\n", jac->smoothsweeps));
1459:       if (jac->smoothtype == 0) {
1460:         PetscStackCallExternalVoid("hypre_ParAMGDataILUType", indx = hypre_ParAMGDataILUType(amg_data));
1461:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU type              %s (%" PetscInt_FMT ")\n", HYPREILUType[indx], indx));
1462:         PetscStackCallExternalVoid("hypre_ParAMGDataILULevel", indx = hypre_ParAMGDataILULevel(amg_data));
1463:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU level             %" PetscInt_FMT "\n", indx));
1464:         PetscStackCallExternalVoid("hypre_ParAMGDataILUMaxIter", indx = hypre_ParAMGDataILUMaxIter(amg_data));
1465:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU max iterations    %" PetscInt_FMT "\n", indx));
1466:         PetscStackCallExternalVoid("hypre_ParAMGDataILUMaxRowNnz", indx = hypre_ParAMGDataILUMaxRowNnz(amg_data));
1467:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU max NNZ per row   %" PetscInt_FMT "\n", indx));
1468:         PetscStackCallExternalVoid("hypre_ParAMGDataILUTriSolve", indx = hypre_ParAMGDataILUTriSolve(amg_data));
1469:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU triangular solve  %" PetscInt_FMT "\n", indx));
1470:         PetscStackCallExternalVoid("hypre_ParAMGDataTol", val = hypre_ParAMGDataTol(amg_data));
1471:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU tolerance         %e\n", val));
1472:         PetscStackCallExternalVoid("hypre_ParAMGDataILUDroptol", val = hypre_ParAMGDataILUDroptol(amg_data));
1473:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU drop tolerance    %e\n", val));
1474:         PetscStackCallExternalVoid("hypre_ParAMGDataILULocalReordering", indx = hypre_ParAMGDataILULocalReordering(amg_data));
1475:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU local reordering  %" PetscInt_FMT "\n", indx));
1476:         PetscStackCallExternalVoid("hypre_ParAMGDataILULowerJacobiIters", indx = hypre_ParAMGDataILULowerJacobiIters(amg_data));
1477:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU lower Jacobi iterations  %" PetscInt_FMT "\n", indx));
1478:         PetscStackCallExternalVoid("hypre_ParAMGDataILUUpperJacobiIters", indx = hypre_ParAMGDataILUUpperJacobiIters(amg_data));
1479:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU upper Jacobi iterations  %" PetscInt_FMT "\n", indx));
1480:         PetscStackCallExternalVoid("hypre_ParAMGDataPrintLevel", indx = hypre_ParAMGDataPrintLevel(amg_data));
1481:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU print level      %" PetscInt_FMT "\n", indx));
1482:         PetscStackCallExternalVoid("hypre_ParAMGDataLogging", indx = hypre_ParAMGDataLogging(amg_data));
1483:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU logging level    %" PetscInt_FMT "\n", indx));
1484:         PetscStackCallExternalVoid("hypre_ParAMGDataILUIterSetupType", indx = hypre_ParAMGDataILUIterSetupType(amg_data));
1485:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU iterative setup type           %s (%" PetscInt_FMT ")\n", HYPREILUIterSetup[indx], indx));
1486:         PetscStackCallExternalVoid("hypre_ParAMGDataILUIterSetupOption", indx = hypre_ParAMGDataILUIterSetupOption(amg_data));
1487:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU iterative setup option         %" PetscInt_FMT "\n", indx));
1488:         PetscStackCallExternalVoid("hypre_ParAMGDataILUIterSetupMaxIter", indx = hypre_ParAMGDataILUIterSetupMaxIter(amg_data));
1489:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU iterative setup max iterations %" PetscInt_FMT "\n", indx));
1490:         PetscStackCallExternalVoid("hypre_ParAMGDataILUIterSetupTolerance", val = hypre_ParAMGDataILUIterSetupTolerance(amg_data));
1491:         PetscCall(PetscViewerASCIIPrintf(viewer, "    ILU iterative setup tolerance      %e\n", val));
1492:       }
1493:     } else {
1494:       PetscCall(PetscViewerASCIIPrintf(viewer, "    Not using more complex smoothers.\n"));
1495:     }
1496:     if (jac->smoothtype == 3) {
1497:       PetscCall(PetscViewerASCIIPrintf(viewer, "    Euclid ILU(k) levels %" PetscInt_FMT "\n", jac->eu_level));
1498:       PetscCall(PetscViewerASCIIPrintf(viewer, "    Euclid ILU(k) drop tolerance %g\n", (double)jac->eu_droptolerance));
1499:       PetscCall(PetscViewerASCIIPrintf(viewer, "    Euclid ILU use Block-Jacobi? %" PetscInt_FMT "\n", jac->eu_bj));
1500:     }
1501:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Measure type        %s\n", HYPREBoomerAMGMeasureType[jac->measuretype]));
1502:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Coarsen type        %s\n", HYPREBoomerAMGCoarsenType[jac->coarsentype]));
1503:     PetscCall(PetscViewerASCIIPrintf(viewer, "    Interpolation type  %s\n", jac->interptype != 100 ? HYPREBoomerAMGInterpType[jac->interptype] : "1pt"));
1504:     if (jac->nodal_coarsening) PetscCall(PetscViewerASCIIPrintf(viewer, "    Using nodal coarsening with HYPRE_BOOMERAMGSetNodal() %" PetscInt_FMT "\n", jac->nodal_coarsening));
1505:     if (jac->vec_interp_variant) {
1506:       PetscCall(PetscViewerASCIIPrintf(viewer, "    HYPRE_BoomerAMGSetInterpVecVariant() %" PetscInt_FMT "\n", jac->vec_interp_variant));
1507:       PetscCall(PetscViewerASCIIPrintf(viewer, "    HYPRE_BoomerAMGSetInterpVecQMax() %" PetscInt_FMT "\n", jac->vec_interp_qmax));
1508:       PetscCall(PetscViewerASCIIPrintf(viewer, "    HYPRE_BoomerAMGSetSmoothInterpVectors() %d\n", jac->vec_interp_smooth));
1509:     }
1510:     if (jac->nodal_relax) PetscCall(PetscViewerASCIIPrintf(viewer, "    Using nodal relaxation via Schwarz smoothing on levels %" PetscInt_FMT "\n", jac->nodal_relax_levels));
1511: #if PETSC_PKG_HYPRE_VERSION_GE(2, 23, 0)
1512:     PetscCall(PetscViewerASCIIPrintf(viewer, "    SpGEMM type         %s\n", jac->spgemm_type));
1513: #else
1514:     PetscCall(PetscViewerASCIIPrintf(viewer, "    SpGEMM type         %s\n", "hypre"));
1515: #endif
1516:     /* AIR */
1517:     if (jac->Rtype) {
1518:       PetscCall(PetscViewerASCIIPrintf(viewer, "    Using approximate ideal restriction type %" PetscInt_FMT "\n", jac->Rtype));
1519:       PetscCall(PetscViewerASCIIPrintf(viewer, "      Threshold for R %g\n", (double)jac->Rstrongthreshold));
1520:       PetscCall(PetscViewerASCIIPrintf(viewer, "      Filter for R %g\n", (double)jac->Rfilterthreshold));
1521:       PetscCall(PetscViewerASCIIPrintf(viewer, "      A drop tolerance %g\n", (double)jac->Adroptol));
1522:       PetscCall(PetscViewerASCIIPrintf(viewer, "      A drop type %" PetscInt_FMT "\n", jac->Adroptype));
1523:     }
1524:   }
1525:   PetscFunctionReturn(PETSC_SUCCESS);
1526: }

1528: static PetscErrorCode PCSetFromOptions_HYPRE_ParaSails(PC pc, PetscOptionItems PetscOptionsObject)
1529: {
1530:   PC_HYPRE   *jac = (PC_HYPRE *)pc->data;
1531:   PetscInt    indx;
1532:   PetscBool   flag;
1533:   const char *symtlist[] = {"nonsymmetric", "SPD", "nonsymmetric,SPD"};

1535:   PetscFunctionBegin;
1536:   PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE ParaSails Options");
1537:   PetscCall(PetscOptionsInt("-pc_hypre_parasails_nlevels", "Number of number of levels", "None", jac->nlevels, &jac->nlevels, 0));
1538:   PetscCall(PetscOptionsReal("-pc_hypre_parasails_thresh", "Threshold", "None", jac->threshold, &jac->threshold, &flag));
1539:   if (flag) PetscCallExternal(HYPRE_ParaSailsSetParams, jac->hsolver, jac->threshold, jac->nlevels);

1541:   PetscCall(PetscOptionsReal("-pc_hypre_parasails_filter", "filter", "None", jac->filter, &jac->filter, &flag));
1542:   if (flag) PetscCallExternal(HYPRE_ParaSailsSetFilter, jac->hsolver, jac->filter);

1544:   PetscCall(PetscOptionsReal("-pc_hypre_parasails_loadbal", "Load balance", "None", jac->loadbal, &jac->loadbal, &flag));
1545:   if (flag) PetscCallExternal(HYPRE_ParaSailsSetLoadbal, jac->hsolver, jac->loadbal);

1547:   PetscCall(PetscOptionsBool("-pc_hypre_parasails_logging", "Print info to screen", "None", (PetscBool)jac->logging, (PetscBool *)&jac->logging, &flag));
1548:   if (flag) PetscCallExternal(HYPRE_ParaSailsSetLogging, jac->hsolver, jac->logging);

1550:   PetscCall(PetscOptionsBool("-pc_hypre_parasails_reuse", "Reuse nonzero pattern in preconditioner", "None", (PetscBool)jac->ruse, (PetscBool *)&jac->ruse, &flag));
1551:   if (flag) PetscCallExternal(HYPRE_ParaSailsSetReuse, jac->hsolver, jac->ruse);

1553:   PetscCall(PetscOptionsEList("-pc_hypre_parasails_sym", "Symmetry of matrix and preconditioner", "None", symtlist, PETSC_STATIC_ARRAY_LENGTH(symtlist), symtlist[0], &indx, &flag));
1554:   if (flag) {
1555:     jac->symt = indx;
1556:     PetscCallExternal(HYPRE_ParaSailsSetSym, jac->hsolver, jac->symt);
1557:   }

1559:   PetscOptionsHeadEnd();
1560:   PetscFunctionReturn(PETSC_SUCCESS);
1561: }

1563: static PetscErrorCode PCView_HYPRE_ParaSails(PC pc, PetscViewer viewer)
1564: {
1565:   PC_HYPRE   *jac = (PC_HYPRE *)pc->data;
1566:   PetscBool   isascii;
1567:   const char *symt = 0;

1569:   PetscFunctionBegin;
1570:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
1571:   if (isascii) {
1572:     PetscCall(PetscViewerASCIIPrintf(viewer, "  HYPRE ParaSails preconditioning\n"));
1573:     PetscCall(PetscViewerASCIIPrintf(viewer, "    nlevels %" PetscInt_FMT "\n", jac->nlevels));
1574:     PetscCall(PetscViewerASCIIPrintf(viewer, "    threshold %g\n", (double)jac->threshold));
1575:     PetscCall(PetscViewerASCIIPrintf(viewer, "    filter %g\n", (double)jac->filter));
1576:     PetscCall(PetscViewerASCIIPrintf(viewer, "    load balance %g\n", (double)jac->loadbal));
1577:     PetscCall(PetscViewerASCIIPrintf(viewer, "    reuse nonzero structure %s\n", PetscBools[jac->ruse]));
1578:     PetscCall(PetscViewerASCIIPrintf(viewer, "    print info to screen %s\n", PetscBools[jac->logging]));
1579:     if (!jac->symt) symt = "nonsymmetric matrix and preconditioner";
1580:     else if (jac->symt == 1) symt = "SPD matrix and preconditioner";
1581:     else if (jac->symt == 2) symt = "nonsymmetric matrix but SPD preconditioner";
1582:     else SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONG, "Unknown HYPRE ParaSails symmetric option %" PetscInt_FMT, jac->symt);
1583:     PetscCall(PetscViewerASCIIPrintf(viewer, "    %s\n", symt));
1584:   }
1585:   PetscFunctionReturn(PETSC_SUCCESS);
1586: }

1588: static PetscErrorCode PCSetFromOptions_HYPRE_AMS(PC pc, PetscOptionItems PetscOptionsObject)
1589: {
1590:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1591:   PetscInt  n;
1592:   PetscBool flag, flag2, flag3, flag4;

1594:   PetscFunctionBegin;
1595:   PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE AMS Options");
1596:   PetscCall(PetscOptionsInt("-pc_hypre_ams_print_level", "Debugging output level for AMS", "None", jac->as_print, &jac->as_print, &flag));
1597:   if (flag) PetscCallExternal(HYPRE_AMSSetPrintLevel, jac->hsolver, jac->as_print);
1598:   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));
1599:   if (flag) PetscCallExternal(HYPRE_AMSSetMaxIter, jac->hsolver, jac->as_max_iter);
1600:   PetscCall(PetscOptionsInt("-pc_hypre_ams_cycle_type", "Cycle type for AMS multigrid", "None", jac->ams_cycle_type, &jac->ams_cycle_type, &flag));
1601:   if (flag) PetscCallExternal(HYPRE_AMSSetCycleType, jac->hsolver, jac->ams_cycle_type);
1602:   PetscCall(PetscOptionsReal("-pc_hypre_ams_tol", "Error tolerance for AMS multigrid", "None", jac->as_tol, &jac->as_tol, &flag));
1603:   if (flag) PetscCallExternal(HYPRE_AMSSetTol, jac->hsolver, jac->as_tol);
1604:   PetscCall(PetscOptionsInt("-pc_hypre_ams_relax_type", "Relaxation type for AMS smoother", "None", jac->as_relax_type, &jac->as_relax_type, &flag));
1605:   PetscCall(PetscOptionsInt("-pc_hypre_ams_relax_times", "Number of relaxation steps for AMS smoother", "None", jac->as_relax_times, &jac->as_relax_times, &flag2));
1606:   PetscCall(PetscOptionsReal("-pc_hypre_ams_relax_weight", "Relaxation weight for AMS smoother", "None", jac->as_relax_weight, &jac->as_relax_weight, &flag3));
1607:   PetscCall(PetscOptionsReal("-pc_hypre_ams_omega", "SSOR coefficient for AMS smoother", "None", jac->as_omega, &jac->as_omega, &flag4));
1608:   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);
1609:   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));
1610:   n = 5;
1611:   PetscCall(PetscOptionsIntArray("-pc_hypre_ams_amg_alpha_options", "AMG options for vector Poisson", "None", jac->as_amg_alpha_opts, &n, &flag2));
1612:   if (flag || flag2) {
1613:     PetscCallExternal(HYPRE_AMSSetAlphaAMGOptions, jac->hsolver, jac->as_amg_alpha_opts[0], /* AMG coarsen type */
1614:                       jac->as_amg_alpha_opts[1],                                            /* AMG agg_levels */
1615:                       jac->as_amg_alpha_opts[2],                                            /* AMG relax_type */
1616:                       jac->as_amg_alpha_theta, jac->as_amg_alpha_opts[3],                   /* AMG interp_type */
1617:                       jac->as_amg_alpha_opts[4]);                                           /* AMG Pmax */
1618:   }
1619:   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));
1620:   n = 5;
1621:   PetscCall(PetscOptionsIntArray("-pc_hypre_ams_amg_beta_options", "AMG options for scalar Poisson solver", "None", jac->as_amg_beta_opts, &n, &flag2));
1622:   if (flag || flag2) {
1623:     PetscCallExternal(HYPRE_AMSSetBetaAMGOptions, jac->hsolver, jac->as_amg_beta_opts[0], /* AMG coarsen type */
1624:                       jac->as_amg_beta_opts[1],                                           /* AMG agg_levels */
1625:                       jac->as_amg_beta_opts[2],                                           /* AMG relax_type */
1626:                       jac->as_amg_beta_theta, jac->as_amg_beta_opts[3],                   /* AMG interp_type */
1627:                       jac->as_amg_beta_opts[4]);                                          /* AMG Pmax */
1628:   }
1629:   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));
1630:   if (flag) { /* override HYPRE's default only if the options is used */
1631:     PetscCallExternal(HYPRE_AMSSetProjectionFrequency, jac->hsolver, jac->ams_proj_freq);
1632:   }
1633:   PetscOptionsHeadEnd();
1634:   PetscFunctionReturn(PETSC_SUCCESS);
1635: }

1637: static PetscErrorCode PCView_HYPRE_AMS(PC pc, PetscViewer viewer)
1638: {
1639:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1640:   PetscBool isascii;

1642:   PetscFunctionBegin;
1643:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
1644:   if (isascii) {
1645:     PetscCall(PetscViewerASCIIPrintf(viewer, "  HYPRE AMS preconditioning\n"));
1646:     PetscCall(PetscViewerASCIIPrintf(viewer, "    subspace iterations per application %" PetscInt_FMT "\n", jac->as_max_iter));
1647:     PetscCall(PetscViewerASCIIPrintf(viewer, "    subspace cycle type %" PetscInt_FMT "\n", jac->ams_cycle_type));
1648:     PetscCall(PetscViewerASCIIPrintf(viewer, "    subspace iteration tolerance %g\n", (double)jac->as_tol));
1649:     PetscCall(PetscViewerASCIIPrintf(viewer, "    smoother type %" PetscInt_FMT "\n", jac->as_relax_type));
1650:     PetscCall(PetscViewerASCIIPrintf(viewer, "    number of smoothing steps %" PetscInt_FMT "\n", jac->as_relax_times));
1651:     PetscCall(PetscViewerASCIIPrintf(viewer, "    smoother weight %g\n", (double)jac->as_relax_weight));
1652:     PetscCall(PetscViewerASCIIPrintf(viewer, "    smoother omega %g\n", (double)jac->as_omega));
1653:     if (jac->alpha_Poisson) {
1654:       PetscCall(PetscViewerASCIIPrintf(viewer, "    vector Poisson solver (passed in by user)\n"));
1655:     } else {
1656:       PetscCall(PetscViewerASCIIPrintf(viewer, "    vector Poisson solver (computed) \n"));
1657:     }
1658:     PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG coarsening type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[0]));
1659:     PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG levels of aggressive coarsening %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[1]));
1660:     PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG relaxation type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[2]));
1661:     PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG interpolation type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[3]));
1662:     PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG max nonzero elements in interpolation rows %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[4]));
1663:     PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG strength threshold %g\n", (double)jac->as_amg_alpha_theta));
1664:     if (!jac->ams_beta_is_zero) {
1665:       if (jac->beta_Poisson) {
1666:         PetscCall(PetscViewerASCIIPrintf(viewer, "    scalar Poisson solver (passed in by user)\n"));
1667:       } else {
1668:         PetscCall(PetscViewerASCIIPrintf(viewer, "    scalar Poisson solver (computed) \n"));
1669:       }
1670:       PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG coarsening type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[0]));
1671:       PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG levels of aggressive coarsening %" PetscInt_FMT "\n", jac->as_amg_beta_opts[1]));
1672:       PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG relaxation type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[2]));
1673:       PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG interpolation type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[3]));
1674:       PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG max nonzero elements in interpolation rows %" PetscInt_FMT "\n", jac->as_amg_beta_opts[4]));
1675:       PetscCall(PetscViewerASCIIPrintf(viewer, "        boomerAMG strength threshold %g\n", (double)jac->as_amg_beta_theta));
1676:       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));
1677:     } else {
1678:       PetscCall(PetscViewerASCIIPrintf(viewer, "    scalar Poisson solver not used (zero-conductivity everywhere) \n"));
1679:     }
1680:   }
1681:   PetscFunctionReturn(PETSC_SUCCESS);
1682: }

1684: static PetscErrorCode PCSetFromOptions_HYPRE_ADS(PC pc, PetscOptionItems PetscOptionsObject)
1685: {
1686:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1687:   PetscInt  n;
1688:   PetscBool flag, flag2, flag3, flag4;

1690:   PetscFunctionBegin;
1691:   PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE ADS Options");
1692:   PetscCall(PetscOptionsInt("-pc_hypre_ads_print_level", "Debugging output level for ADS", "None", jac->as_print, &jac->as_print, &flag));
1693:   if (flag) PetscCallExternal(HYPRE_ADSSetPrintLevel, jac->hsolver, jac->as_print);
1694:   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));
1695:   if (flag) PetscCallExternal(HYPRE_ADSSetMaxIter, jac->hsolver, jac->as_max_iter);
1696:   PetscCall(PetscOptionsInt("-pc_hypre_ads_cycle_type", "Cycle type for ADS multigrid", "None", jac->ads_cycle_type, &jac->ads_cycle_type, &flag));
1697:   if (flag) PetscCallExternal(HYPRE_ADSSetCycleType, jac->hsolver, jac->ads_cycle_type);
1698:   PetscCall(PetscOptionsReal("-pc_hypre_ads_tol", "Error tolerance for ADS multigrid", "None", jac->as_tol, &jac->as_tol, &flag));
1699:   if (flag) PetscCallExternal(HYPRE_ADSSetTol, jac->hsolver, jac->as_tol);
1700:   PetscCall(PetscOptionsInt("-pc_hypre_ads_relax_type", "Relaxation type for ADS smoother", "None", jac->as_relax_type, &jac->as_relax_type, &flag));
1701:   PetscCall(PetscOptionsInt("-pc_hypre_ads_relax_times", "Number of relaxation steps for ADS smoother", "None", jac->as_relax_times, &jac->as_relax_times, &flag2));
1702:   PetscCall(PetscOptionsReal("-pc_hypre_ads_relax_weight", "Relaxation weight for ADS smoother", "None", jac->as_relax_weight, &jac->as_relax_weight, &flag3));
1703:   PetscCall(PetscOptionsReal("-pc_hypre_ads_omega", "SSOR coefficient for ADS smoother", "None", jac->as_omega, &jac->as_omega, &flag4));
1704:   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);
1705:   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));
1706:   n = 5;
1707:   PetscCall(PetscOptionsIntArray("-pc_hypre_ads_ams_options", "AMG options for AMS solver inside ADS", "None", jac->as_amg_alpha_opts, &n, &flag2));
1708:   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));
1709:   if (flag || flag2 || flag3) {
1710:     PetscCallExternal(HYPRE_ADSSetAMSOptions, jac->hsolver, jac->ams_cycle_type, /* AMS cycle type */
1711:                       jac->as_amg_alpha_opts[0],                                 /* AMG coarsen type */
1712:                       jac->as_amg_alpha_opts[1],                                 /* AMG agg_levels */
1713:                       jac->as_amg_alpha_opts[2],                                 /* AMG relax_type */
1714:                       jac->as_amg_alpha_theta, jac->as_amg_alpha_opts[3],        /* AMG interp_type */
1715:                       jac->as_amg_alpha_opts[4]);                                /* AMG Pmax */
1716:   }
1717:   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));
1718:   n = 5;
1719:   PetscCall(PetscOptionsIntArray("-pc_hypre_ads_amg_options", "AMG options for vector AMG solver inside ADS", "None", jac->as_amg_beta_opts, &n, &flag2));
1720:   if (flag || flag2) {
1721:     PetscCallExternal(HYPRE_ADSSetAMGOptions, jac->hsolver, jac->as_amg_beta_opts[0], /* AMG coarsen type */
1722:                       jac->as_amg_beta_opts[1],                                       /* AMG agg_levels */
1723:                       jac->as_amg_beta_opts[2],                                       /* AMG relax_type */
1724:                       jac->as_amg_beta_theta, jac->as_amg_beta_opts[3],               /* AMG interp_type */
1725:                       jac->as_amg_beta_opts[4]);                                      /* AMG Pmax */
1726:   }
1727:   PetscOptionsHeadEnd();
1728:   PetscFunctionReturn(PETSC_SUCCESS);
1729: }

1731: static PetscErrorCode PCView_HYPRE_ADS(PC pc, PetscViewer viewer)
1732: {
1733:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1734:   PetscBool isascii;

1736:   PetscFunctionBegin;
1737:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
1738:   if (isascii) {
1739:     PetscCall(PetscViewerASCIIPrintf(viewer, "  HYPRE ADS preconditioning\n"));
1740:     PetscCall(PetscViewerASCIIPrintf(viewer, "    subspace iterations per application %" PetscInt_FMT "\n", jac->as_max_iter));
1741:     PetscCall(PetscViewerASCIIPrintf(viewer, "    subspace cycle type %" PetscInt_FMT "\n", jac->ads_cycle_type));
1742:     PetscCall(PetscViewerASCIIPrintf(viewer, "    subspace iteration tolerance %g\n", (double)jac->as_tol));
1743:     PetscCall(PetscViewerASCIIPrintf(viewer, "    smoother type %" PetscInt_FMT "\n", jac->as_relax_type));
1744:     PetscCall(PetscViewerASCIIPrintf(viewer, "    number of smoothing steps %" PetscInt_FMT "\n", jac->as_relax_times));
1745:     PetscCall(PetscViewerASCIIPrintf(viewer, "    smoother weight %g\n", (double)jac->as_relax_weight));
1746:     PetscCall(PetscViewerASCIIPrintf(viewer, "    smoother omega %g\n", (double)jac->as_omega));
1747:     PetscCall(PetscViewerASCIIPrintf(viewer, "    AMS solver using boomerAMG\n"));
1748:     PetscCall(PetscViewerASCIIPrintf(viewer, "        subspace cycle type %" PetscInt_FMT "\n", jac->ams_cycle_type));
1749:     PetscCall(PetscViewerASCIIPrintf(viewer, "        coarsening type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[0]));
1750:     PetscCall(PetscViewerASCIIPrintf(viewer, "        levels of aggressive coarsening %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[1]));
1751:     PetscCall(PetscViewerASCIIPrintf(viewer, "        relaxation type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[2]));
1752:     PetscCall(PetscViewerASCIIPrintf(viewer, "        interpolation type %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[3]));
1753:     PetscCall(PetscViewerASCIIPrintf(viewer, "        max nonzero elements in interpolation rows %" PetscInt_FMT "\n", jac->as_amg_alpha_opts[4]));
1754:     PetscCall(PetscViewerASCIIPrintf(viewer, "        strength threshold %g\n", (double)jac->as_amg_alpha_theta));
1755:     PetscCall(PetscViewerASCIIPrintf(viewer, "    vector Poisson solver using boomerAMG\n"));
1756:     PetscCall(PetscViewerASCIIPrintf(viewer, "        coarsening type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[0]));
1757:     PetscCall(PetscViewerASCIIPrintf(viewer, "        levels of aggressive coarsening %" PetscInt_FMT "\n", jac->as_amg_beta_opts[1]));
1758:     PetscCall(PetscViewerASCIIPrintf(viewer, "        relaxation type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[2]));
1759:     PetscCall(PetscViewerASCIIPrintf(viewer, "        interpolation type %" PetscInt_FMT "\n", jac->as_amg_beta_opts[3]));
1760:     PetscCall(PetscViewerASCIIPrintf(viewer, "        max nonzero elements in interpolation rows %" PetscInt_FMT "\n", jac->as_amg_beta_opts[4]));
1761:     PetscCall(PetscViewerASCIIPrintf(viewer, "        strength threshold %g\n", (double)jac->as_amg_beta_theta));
1762:   }
1763:   PetscFunctionReturn(PETSC_SUCCESS);
1764: }

1766: static PetscErrorCode PCHYPRESetDiscreteGradient_HYPRE(PC pc, Mat G)
1767: {
1768:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1769:   PetscBool ishypre;

1771:   PetscFunctionBegin;
1772:   PetscCall(PetscObjectTypeCompare((PetscObject)G, MATHYPRE, &ishypre));
1773:   if (ishypre) {
1774:     PetscCall(PetscObjectReference((PetscObject)G));
1775:     PetscCall(MatDestroy(&jac->G));
1776:     jac->G = G;
1777:   } else {
1778:     PetscCall(MatDestroy(&jac->G));
1779:     PetscCall(MatConvert(G, MATHYPRE, MAT_INITIAL_MATRIX, &jac->G));
1780:   }
1781:   PetscFunctionReturn(PETSC_SUCCESS);
1782: }

1784: /*@
1785:   PCHYPRESetDiscreteGradient - Set discrete gradient matrix for `PCHYPRE` type of ams or ads

1787:   Collective

1789:   Input Parameters:
1790: + pc - the preconditioning context
1791: - G  - the discrete gradient

1793:   Level: intermediate

1795:   Notes:
1796:   G should have as many rows as the number of edges and as many columns as the number of vertices in the mesh

1798:   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

1800:   Developer Notes:
1801:   This automatically converts the matrix to `MATHYPRE` if it is not already of that type

1803: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteCurl()`
1804: @*/
1805: PetscErrorCode PCHYPRESetDiscreteGradient(PC pc, Mat G)
1806: {
1807:   PetscFunctionBegin;
1810:   PetscCheckSameComm(pc, 1, G, 2);
1811:   PetscTryMethod(pc, "PCHYPRESetDiscreteGradient_C", (PC, Mat), (pc, G));
1812:   PetscFunctionReturn(PETSC_SUCCESS);
1813: }

1815: static PetscErrorCode PCHYPRESetDiscreteCurl_HYPRE(PC pc, Mat C)
1816: {
1817:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1818:   PetscBool ishypre;

1820:   PetscFunctionBegin;
1821:   PetscCall(PetscObjectTypeCompare((PetscObject)C, MATHYPRE, &ishypre));
1822:   if (ishypre) {
1823:     PetscCall(PetscObjectReference((PetscObject)C));
1824:     PetscCall(MatDestroy(&jac->C));
1825:     jac->C = C;
1826:   } else {
1827:     PetscCall(MatDestroy(&jac->C));
1828:     PetscCall(MatConvert(C, MATHYPRE, MAT_INITIAL_MATRIX, &jac->C));
1829:   }
1830:   PetscFunctionReturn(PETSC_SUCCESS);
1831: }

1833: /*@
1834:   PCHYPRESetDiscreteCurl - Set discrete curl matrix for `PCHYPRE` type of ads

1836:   Collective

1838:   Input Parameters:
1839: + pc - the preconditioning context
1840: - C  - the discrete curl

1842:   Level: intermediate

1844:   Notes:
1845:   C should have as many rows as the number of faces and as many columns as the number of edges in the mesh

1847:   Each row of C 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

1849:   Developer Notes:
1850:   This automatically converts the matrix to `MATHYPRE` if it is not already of that type

1852:   If this is only for  `PCHYPRE` type of ads it should be called `PCHYPREADSSetDiscreteCurl()`

1854: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteGradient()`
1855: @*/
1856: PetscErrorCode PCHYPRESetDiscreteCurl(PC pc, Mat C)
1857: {
1858:   PetscFunctionBegin;
1861:   PetscCheckSameComm(pc, 1, C, 2);
1862:   PetscTryMethod(pc, "PCHYPRESetDiscreteCurl_C", (PC, Mat), (pc, C));
1863:   PetscFunctionReturn(PETSC_SUCCESS);
1864: }

1866: static PetscErrorCode PCHYPRESetInterpolations_HYPRE(PC pc, PetscInt dim, Mat RT_PiFull, Mat RT_Pi[], Mat ND_PiFull, Mat ND_Pi[])
1867: {
1868:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1869:   PetscBool ishypre;
1870:   PetscInt  i;

1872:   PetscFunctionBegin;
1873:   PetscCall(MatDestroy(&jac->RT_PiFull));
1874:   PetscCall(MatDestroy(&jac->ND_PiFull));
1875:   for (i = 0; i < 3; ++i) {
1876:     PetscCall(MatDestroy(&jac->RT_Pi[i]));
1877:     PetscCall(MatDestroy(&jac->ND_Pi[i]));
1878:   }

1880:   jac->dim = dim;
1881:   if (RT_PiFull) {
1882:     PetscCall(PetscObjectTypeCompare((PetscObject)RT_PiFull, MATHYPRE, &ishypre));
1883:     if (ishypre) {
1884:       PetscCall(PetscObjectReference((PetscObject)RT_PiFull));
1885:       jac->RT_PiFull = RT_PiFull;
1886:     } else {
1887:       PetscCall(MatConvert(RT_PiFull, MATHYPRE, MAT_INITIAL_MATRIX, &jac->RT_PiFull));
1888:     }
1889:   }
1890:   if (RT_Pi) {
1891:     for (i = 0; i < dim; ++i) {
1892:       if (RT_Pi[i]) {
1893:         PetscCall(PetscObjectTypeCompare((PetscObject)RT_Pi[i], MATHYPRE, &ishypre));
1894:         if (ishypre) {
1895:           PetscCall(PetscObjectReference((PetscObject)RT_Pi[i]));
1896:           jac->RT_Pi[i] = RT_Pi[i];
1897:         } else {
1898:           PetscCall(MatConvert(RT_Pi[i], MATHYPRE, MAT_INITIAL_MATRIX, &jac->RT_Pi[i]));
1899:         }
1900:       }
1901:     }
1902:   }
1903:   if (ND_PiFull) {
1904:     PetscCall(PetscObjectTypeCompare((PetscObject)ND_PiFull, MATHYPRE, &ishypre));
1905:     if (ishypre) {
1906:       PetscCall(PetscObjectReference((PetscObject)ND_PiFull));
1907:       jac->ND_PiFull = ND_PiFull;
1908:     } else {
1909:       PetscCall(MatConvert(ND_PiFull, MATHYPRE, MAT_INITIAL_MATRIX, &jac->ND_PiFull));
1910:     }
1911:   }
1912:   if (ND_Pi) {
1913:     for (i = 0; i < dim; ++i) {
1914:       if (ND_Pi[i]) {
1915:         PetscCall(PetscObjectTypeCompare((PetscObject)ND_Pi[i], MATHYPRE, &ishypre));
1916:         if (ishypre) {
1917:           PetscCall(PetscObjectReference((PetscObject)ND_Pi[i]));
1918:           jac->ND_Pi[i] = ND_Pi[i];
1919:         } else {
1920:           PetscCall(MatConvert(ND_Pi[i], MATHYPRE, MAT_INITIAL_MATRIX, &jac->ND_Pi[i]));
1921:         }
1922:       }
1923:     }
1924:   }
1925:   PetscFunctionReturn(PETSC_SUCCESS);
1926: }

1928: /*@
1929:   PCHYPRESetInterpolations - Set interpolation matrices for `PCHYPRE` type of ams or ads

1931:   Collective

1933:   Input Parameters:
1934: + pc        - the preconditioning context
1935: . dim       - the dimension of the problem, only used in AMS
1936: . RT_PiFull - Raviart-Thomas interpolation matrix
1937: . RT_Pi     - x/y/z component of Raviart-Thomas interpolation matrix
1938: . ND_PiFull - Nedelec interpolation matrix
1939: - ND_Pi     - x/y/z component of Nedelec interpolation matrix

1941:   Level: intermediate

1943:   Notes:
1944:   For AMS, only Nedelec interpolation matrices are needed, the Raviart-Thomas interpolation matrices can be set to NULL.

1946:   For ADS, both type of interpolation matrices are needed.

1948:   Developer Notes:
1949:   This automatically converts the matrix to `MATHYPRE` if it is not already of that type

1951: .seealso: [](ch_ksp), `PCHYPRE`
1952: @*/
1953: PetscErrorCode PCHYPRESetInterpolations(PC pc, PetscInt dim, Mat RT_PiFull, Mat RT_Pi[], Mat ND_PiFull, Mat ND_Pi[])
1954: {
1955:   PetscInt i;

1957:   PetscFunctionBegin;
1959:   if (RT_PiFull) {
1961:     PetscCheckSameComm(pc, 1, RT_PiFull, 3);
1962:   }
1963:   if (RT_Pi) {
1964:     PetscAssertPointer(RT_Pi, 4);
1965:     for (i = 0; i < dim; ++i) {
1966:       if (RT_Pi[i]) {
1968:         PetscCheckSameComm(pc, 1, RT_Pi[i], 4);
1969:       }
1970:     }
1971:   }
1972:   if (ND_PiFull) {
1974:     PetscCheckSameComm(pc, 1, ND_PiFull, 5);
1975:   }
1976:   if (ND_Pi) {
1977:     PetscAssertPointer(ND_Pi, 6);
1978:     for (i = 0; i < dim; ++i) {
1979:       if (ND_Pi[i]) {
1981:         PetscCheckSameComm(pc, 1, ND_Pi[i], 6);
1982:       }
1983:     }
1984:   }
1985:   PetscTryMethod(pc, "PCHYPRESetInterpolations_C", (PC, PetscInt, Mat, Mat[], Mat, Mat[]), (pc, dim, RT_PiFull, RT_Pi, ND_PiFull, ND_Pi));
1986:   PetscFunctionReturn(PETSC_SUCCESS);
1987: }

1989: static PetscErrorCode PCHYPRESetPoissonMatrix_HYPRE(PC pc, Mat A, PetscBool isalpha)
1990: {
1991:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
1992:   PetscBool ishypre;

1994:   PetscFunctionBegin;
1995:   PetscCall(PetscObjectTypeCompare((PetscObject)A, MATHYPRE, &ishypre));
1996:   if (ishypre) {
1997:     if (isalpha) {
1998:       PetscCall(PetscObjectReference((PetscObject)A));
1999:       PetscCall(MatDestroy(&jac->alpha_Poisson));
2000:       jac->alpha_Poisson = A;
2001:     } else {
2002:       if (A) {
2003:         PetscCall(PetscObjectReference((PetscObject)A));
2004:       } else {
2005:         jac->ams_beta_is_zero = PETSC_TRUE;
2006:       }
2007:       PetscCall(MatDestroy(&jac->beta_Poisson));
2008:       jac->beta_Poisson = A;
2009:     }
2010:   } else {
2011:     if (isalpha) {
2012:       PetscCall(MatDestroy(&jac->alpha_Poisson));
2013:       PetscCall(MatConvert(A, MATHYPRE, MAT_INITIAL_MATRIX, &jac->alpha_Poisson));
2014:     } else {
2015:       if (A) {
2016:         PetscCall(MatDestroy(&jac->beta_Poisson));
2017:         PetscCall(MatConvert(A, MATHYPRE, MAT_INITIAL_MATRIX, &jac->beta_Poisson));
2018:       } else {
2019:         PetscCall(MatDestroy(&jac->beta_Poisson));
2020:         jac->ams_beta_is_zero = PETSC_TRUE;
2021:       }
2022:     }
2023:   }
2024:   PetscFunctionReturn(PETSC_SUCCESS);
2025: }

2027: /*@
2028:   PCHYPRESetAlphaPoissonMatrix - Set vector Poisson matrix for `PCHYPRE` of type ams

2030:   Collective

2032:   Input Parameters:
2033: + pc - the preconditioning context
2034: - A  - the matrix

2036:   Level: intermediate

2038:   Note:
2039:   A should be obtained by discretizing the vector valued Poisson problem with linear finite elements

2041:   Developer Notes:
2042:   This automatically converts the matrix to `MATHYPRE` if it is not already of that type

2044:   If this is only for  `PCHYPRE` type of ams it should be called `PCHYPREAMSSetAlphaPoissonMatrix()`

2046: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetBetaPoissonMatrix()`
2047: @*/
2048: PetscErrorCode PCHYPRESetAlphaPoissonMatrix(PC pc, Mat A)
2049: {
2050:   PetscFunctionBegin;
2053:   PetscCheckSameComm(pc, 1, A, 2);
2054:   PetscTryMethod(pc, "PCHYPRESetPoissonMatrix_C", (PC, Mat, PetscBool), (pc, A, PETSC_TRUE));
2055:   PetscFunctionReturn(PETSC_SUCCESS);
2056: }

2058: /*@
2059:   PCHYPRESetBetaPoissonMatrix - Set Poisson matrix for `PCHYPRE` of type ams

2061:   Collective

2063:   Input Parameters:
2064: + pc - the preconditioning context
2065: - A  - the matrix, or NULL to turn it off

2067:   Level: intermediate

2069:   Note:
2070:   A should be obtained by discretizing the Poisson problem with linear finite elements.

2072:   Developer Notes:
2073:   This automatically converts the matrix to `MATHYPRE` if it is not already of that type

2075:   If this is only for  `PCHYPRE` type of ams it should be called `PCHYPREAMSPCHYPRESetBetaPoissonMatrix()`

2077: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetAlphaPoissonMatrix()`
2078: @*/
2079: PetscErrorCode PCHYPRESetBetaPoissonMatrix(PC pc, Mat A)
2080: {
2081:   PetscFunctionBegin;
2083:   if (A) {
2085:     PetscCheckSameComm(pc, 1, A, 2);
2086:   }
2087:   PetscTryMethod(pc, "PCHYPRESetPoissonMatrix_C", (PC, Mat, PetscBool), (pc, A, PETSC_FALSE));
2088:   PetscFunctionReturn(PETSC_SUCCESS);
2089: }

2091: static PetscErrorCode PCHYPRESetEdgeConstantVectors_HYPRE(PC pc, Vec ozz, Vec zoz, Vec zzo)
2092: {
2093:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;

2095:   PetscFunctionBegin;
2096:   /* throw away any vector if already set */
2097:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[0]));
2098:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[1]));
2099:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->constants[2]));
2100:   PetscCall(VecHYPRE_IJVectorCreate(ozz->map, &jac->constants[0]));
2101:   PetscCall(VecHYPRE_IJVectorCopy(ozz, jac->constants[0]));
2102:   PetscCall(VecHYPRE_IJVectorCreate(zoz->map, &jac->constants[1]));
2103:   PetscCall(VecHYPRE_IJVectorCopy(zoz, jac->constants[1]));
2104:   jac->dim = 2;
2105:   if (zzo) {
2106:     PetscCall(VecHYPRE_IJVectorCreate(zzo->map, &jac->constants[2]));
2107:     PetscCall(VecHYPRE_IJVectorCopy(zzo, jac->constants[2]));
2108:     jac->dim++;
2109:   }
2110:   PetscFunctionReturn(PETSC_SUCCESS);
2111: }

2113: /*@
2114:   PCHYPRESetEdgeConstantVectors - Set the representation of the constant vector fields in the edge element basis for `PCHYPRE` of type ams

2116:   Collective

2118:   Input Parameters:
2119: + pc  - the preconditioning context
2120: . ozz - vector representing (1,0,0) (or (1,0) in 2D)
2121: . zoz - vector representing (0,1,0) (or (0,1) in 2D)
2122: - zzo - vector representing (0,0,1) (use NULL in 2D)

2124:   Level: intermediate

2126:   Developer Notes:
2127:   If this is only for  `PCHYPRE` type of ams it should be called `PCHYPREAMSSetEdgeConstantVectors()`

2129: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetAlphaPoissonMatrix()`
2130: @*/
2131: PetscErrorCode PCHYPRESetEdgeConstantVectors(PC pc, Vec ozz, Vec zoz, Vec zzo)
2132: {
2133:   PetscFunctionBegin;
2138:   PetscCheckSameComm(pc, 1, ozz, 2);
2139:   PetscCheckSameComm(pc, 1, zoz, 3);
2140:   if (zzo) PetscCheckSameComm(pc, 1, zzo, 4);
2141:   PetscTryMethod(pc, "PCHYPRESetEdgeConstantVectors_C", (PC, Vec, Vec, Vec), (pc, ozz, zoz, zzo));
2142:   PetscFunctionReturn(PETSC_SUCCESS);
2143: }

2145: static PetscErrorCode PCHYPREAMSSetInteriorNodes_HYPRE(PC pc, Vec interior)
2146: {
2147:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;

2149:   PetscFunctionBegin;
2150:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->interior));
2151:   PetscCall(VecHYPRE_IJVectorCreate(interior->map, &jac->interior));
2152:   PetscCall(VecHYPRE_IJVectorCopy(interior, jac->interior));
2153:   jac->ams_beta_is_zero_part = PETSC_TRUE;
2154:   PetscFunctionReturn(PETSC_SUCCESS);
2155: }

2157: /*@
2158:   PCHYPREAMSSetInteriorNodes - Set the list of interior nodes to a zero-conductivity region for `PCHYPRE` of type ams

2160:   Collective

2162:   Input Parameters:
2163: + pc       - the preconditioning context
2164: - interior - vector. node is interior if its entry in the array is 1.0.

2166:   Level: intermediate

2168:   Note:
2169:   This calls `HYPRE_AMSSetInteriorNodes()`

2171: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetAlphaPoissonMatrix()`
2172: @*/
2173: PetscErrorCode PCHYPREAMSSetInteriorNodes(PC pc, Vec interior)
2174: {
2175:   PetscFunctionBegin;
2178:   PetscCheckSameComm(pc, 1, interior, 2);
2179:   PetscTryMethod(pc, "PCHYPREAMSSetInteriorNodes_C", (PC, Vec), (pc, interior));
2180:   PetscFunctionReturn(PETSC_SUCCESS);
2181: }

2183: static PetscErrorCode PCSetCoordinates_HYPRE(PC pc, PetscInt dim, PetscInt nloc, PetscReal *coords)
2184: {
2185:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
2186:   Vec       tv;
2187:   PetscInt  i;

2189:   PetscFunctionBegin;
2190:   /* throw away any coordinate vector if already set */
2191:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[0]));
2192:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[1]));
2193:   PetscCall(VecHYPRE_IJVectorDestroy(&jac->coords[2]));
2194:   jac->dim = dim;

2196:   /* compute IJ vector for coordinates */
2197:   PetscCall(VecCreate(PetscObjectComm((PetscObject)pc), &tv));
2198:   PetscCall(VecSetType(tv, VECSTANDARD));
2199:   PetscCall(VecSetSizes(tv, nloc, PETSC_DECIDE));
2200:   for (i = 0; i < dim; i++) {
2201:     PetscScalar *array;
2202:     PetscInt     j;

2204:     PetscCall(VecHYPRE_IJVectorCreate(tv->map, &jac->coords[i]));
2205:     PetscCall(VecGetArrayWrite(tv, &array));
2206:     for (j = 0; j < nloc; j++) array[j] = coords[j * dim + i];
2207:     PetscCall(VecRestoreArrayWrite(tv, &array));
2208:     PetscCall(VecHYPRE_IJVectorCopy(tv, jac->coords[i]));
2209:   }
2210:   PetscCall(VecDestroy(&tv));
2211:   PetscFunctionReturn(PETSC_SUCCESS);
2212: }

2214: static PetscErrorCode PCHYPREGetType_HYPRE(PC pc, const char *name[])
2215: {
2216:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;

2218:   PetscFunctionBegin;
2219:   *name = jac->hypre_type;
2220:   PetscFunctionReturn(PETSC_SUCCESS);
2221: }

2223: static PetscErrorCode PCHYPRESetType_HYPRE(PC pc, const char name[])
2224: {
2225:   PC_HYPRE *jac = (PC_HYPRE *)pc->data;
2226:   PetscBool flag;

2228:   PetscFunctionBegin;
2229:   if (jac->hypre_type) {
2230:     PetscCall(PetscStrcmp(jac->hypre_type, name, &flag));
2231:     if (flag) PetscFunctionReturn(PETSC_SUCCESS);
2232:   }

2234:   PetscCall(PCReset_HYPRE(pc));
2235:   PetscCall(PetscFree(jac->hypre_type));
2236:   PetscCall(PetscStrallocpy(name, &jac->hypre_type));

2238:   jac->maxiter         = PETSC_DEFAULT;
2239:   jac->tol             = PETSC_DEFAULT;
2240:   jac->printstatistics = PetscLogPrintInfo;

2242:   PetscCall(PetscStrcmp("ilu", jac->hypre_type, &flag));
2243:   if (flag) {
2244:     PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &jac->comm_hypre));
2245:     PetscCallExternal(HYPRE_ILUCreate, &jac->hsolver);
2246:     pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ILU;
2247:     pc->ops->view           = PCView_HYPRE_ILU;
2248:     jac->destroy            = HYPRE_ILUDestroy;
2249:     jac->setup              = HYPRE_ILUSetup;
2250:     jac->solve              = HYPRE_ILUSolve;
2251:     jac->factorrowsize      = PETSC_DEFAULT;
2252:     PetscFunctionReturn(PETSC_SUCCESS);
2253:   }

2255:   PetscCall(PetscStrcmp("pilut", jac->hypre_type, &flag));
2256:   if (flag) {
2257:     PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &jac->comm_hypre));
2258:     PetscCallExternal(HYPRE_ParCSRPilutCreate, jac->comm_hypre, &jac->hsolver);
2259:     pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Pilut;
2260:     pc->ops->view           = PCView_HYPRE_Pilut;
2261:     jac->destroy            = HYPRE_ParCSRPilutDestroy;
2262:     jac->setup              = HYPRE_ParCSRPilutSetup;
2263:     jac->solve              = HYPRE_ParCSRPilutSolve;
2264:     jac->factorrowsize      = PETSC_DEFAULT;
2265:     PetscFunctionReturn(PETSC_SUCCESS);
2266:   }
2267:   PetscCall(PetscStrcmp("euclid", jac->hypre_type, &flag));
2268:   if (flag) {
2269: #if defined(PETSC_USE_64BIT_INDICES)
2270:     SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_SUP, "Hypre Euclid does not support 64-bit indices");
2271: #endif
2272:     PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &jac->comm_hypre));
2273:     PetscCallExternal(HYPRE_EuclidCreate, jac->comm_hypre, &jac->hsolver);
2274:     pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Euclid;
2275:     pc->ops->view           = PCView_HYPRE_Euclid;
2276:     jac->destroy            = HYPRE_EuclidDestroy;
2277:     jac->setup              = HYPRE_EuclidSetup;
2278:     jac->solve              = HYPRE_EuclidSolve;
2279:     jac->factorrowsize      = PETSC_DEFAULT;
2280:     jac->eu_level           = PETSC_DEFAULT; /* default */
2281:     PetscFunctionReturn(PETSC_SUCCESS);
2282:   }
2283:   PetscCall(PetscStrcmp("parasails", jac->hypre_type, &flag));
2284:   if (flag) {
2285:     PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &jac->comm_hypre));
2286:     PetscCallExternal(HYPRE_ParaSailsCreate, jac->comm_hypre, &jac->hsolver);
2287:     pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ParaSails;
2288:     pc->ops->view           = PCView_HYPRE_ParaSails;
2289:     jac->destroy            = HYPRE_ParaSailsDestroy;
2290:     jac->setup              = HYPRE_ParaSailsSetup;
2291:     jac->solve              = HYPRE_ParaSailsSolve;
2292:     /* initialize */
2293:     jac->nlevels   = 1;
2294:     jac->threshold = .1;
2295:     jac->filter    = .1;
2296:     jac->loadbal   = 0;
2297:     if (PetscLogPrintInfo) jac->logging = (int)PETSC_TRUE;
2298:     else jac->logging = (int)PETSC_FALSE;

2300:     jac->ruse = (int)PETSC_FALSE;
2301:     jac->symt = 0;
2302:     PetscCallExternal(HYPRE_ParaSailsSetParams, jac->hsolver, jac->threshold, jac->nlevels);
2303:     PetscCallExternal(HYPRE_ParaSailsSetFilter, jac->hsolver, jac->filter);
2304:     PetscCallExternal(HYPRE_ParaSailsSetLoadbal, jac->hsolver, jac->loadbal);
2305:     PetscCallExternal(HYPRE_ParaSailsSetLogging, jac->hsolver, jac->logging);
2306:     PetscCallExternal(HYPRE_ParaSailsSetReuse, jac->hsolver, jac->ruse);
2307:     PetscCallExternal(HYPRE_ParaSailsSetSym, jac->hsolver, jac->symt);
2308:     PetscFunctionReturn(PETSC_SUCCESS);
2309:   }
2310:   PetscCall(PetscStrcmp("boomeramg", jac->hypre_type, &flag));
2311:   if (flag) {
2312:     PetscCallExternal(HYPRE_BoomerAMGCreate, &jac->hsolver);
2313:     pc->ops->setfromoptions  = PCSetFromOptions_HYPRE_BoomerAMG;
2314:     pc->ops->view            = PCView_HYPRE_BoomerAMG;
2315:     pc->ops->applytranspose  = PCApplyTranspose_HYPRE_BoomerAMG;
2316:     pc->ops->applyrichardson = PCApplyRichardson_HYPRE_BoomerAMG;
2317:     pc->ops->matapply        = PCMatApply_HYPRE_BoomerAMG;
2318:     PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCGetInterpolations_C", PCGetInterpolations_BoomerAMG));
2319:     PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCGetCoarseOperators_C", PCGetCoarseOperators_BoomerAMG));
2320:     PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPREGetCFMarkers_C", PCHYPREGetCFMarkers_BoomerAMG));
2321:     jac->destroy         = HYPRE_BoomerAMGDestroy;
2322:     jac->setup           = HYPRE_BoomerAMGSetup;
2323:     jac->solve           = HYPRE_BoomerAMGSolve;
2324:     jac->applyrichardson = PETSC_FALSE;
2325:     /* these defaults match the hypre defaults */
2326:     jac->cycletype       = 1;
2327:     jac->maxlevels       = 25;
2328:     jac->maxiter         = 1;
2329:     jac->tol             = 0.0; /* tolerance of zero indicates use as preconditioner (suppresses convergence errors) */
2330:     jac->truncfactor     = 0.0;
2331:     jac->strongthreshold = .25;
2332:     jac->maxrowsum       = .9;
2333:     jac->measuretype     = 0;
2334:     jac->gridsweeps[0] = jac->gridsweeps[1] = jac->gridsweeps[2] = 1;
2335:     jac->smoothtype                                              = -1; /* Not set by default */
2336:     jac->smoothnumlevels                                         = 25;
2337:     jac->eu_level                                                = 0;
2338:     jac->eu_droptolerance                                        = 0;
2339:     jac->eu_bj                                                   = 0;
2340:     jac->relaxweight                                             = 1.0;
2341:     jac->outerrelaxweight                                        = 1.0;
2342:     jac->Rtype                                                   = 0;
2343:     jac->Rstrongthreshold                                        = 0.25;
2344:     jac->Rfilterthreshold                                        = 0.0;
2345:     jac->Adroptype                                               = -1;
2346:     jac->Adroptol                                                = 0.0;
2347:     jac->agg_nl                                                  = 0;
2348:     jac->pmax                                                    = 0;
2349:     jac->truncfactor                                             = 0.0;
2350:     jac->agg_num_paths                                           = 1;
2351:     jac->maxc                                                    = 9;
2352:     jac->minc                                                    = 1;
2353:     jac->nodal_coarsening                                        = 0;
2354:     jac->nodal_coarsening_diag                                   = 0;
2355:     jac->vec_interp_variant                                      = 0;
2356:     jac->vec_interp_qmax                                         = 0;
2357:     jac->vec_interp_smooth                                       = PETSC_FALSE;
2358:     jac->interp_refine                                           = 0;
2359:     jac->nodal_relax                                             = PETSC_FALSE;
2360:     jac->nodal_relax_levels                                      = 1;
2361:     jac->rap2                                                    = 0;
2362:     PetscObjectParameterSetDefault(jac, relaxorder, -1); /* Initialize with invalid value so we can recognize user input */
2363:     PetscFunctionReturn(PETSC_SUCCESS);
2364:   }
2365:   PetscCall(PetscStrcmp("ams", jac->hypre_type, &flag));
2366:   if (flag) {
2367:     PetscCallExternal(HYPRE_AMSCreate, &jac->hsolver);
2368:     pc->ops->setfromoptions = PCSetFromOptions_HYPRE_AMS;
2369:     pc->ops->view           = PCView_HYPRE_AMS;
2370:     jac->destroy            = HYPRE_AMSDestroy;
2371:     jac->setup              = HYPRE_AMSSetup;
2372:     jac->solve              = HYPRE_AMSSolve;
2373:     jac->coords[0]          = NULL;
2374:     jac->coords[1]          = NULL;
2375:     jac->coords[2]          = NULL;
2376:     jac->interior           = NULL;
2377:     /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE AMS */
2378:     jac->as_print       = 0;
2379:     jac->as_max_iter    = 1;  /* used as a preconditioner */
2380:     jac->as_tol         = 0.; /* used as a preconditioner */
2381:     jac->ams_cycle_type = 13;
2382:     /* Smoothing options */
2383:     jac->as_relax_type   = 2;
2384:     jac->as_relax_times  = 1;
2385:     jac->as_relax_weight = 1.0;
2386:     jac->as_omega        = 1.0;
2387:     /* Vector valued Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
2388:     jac->as_amg_alpha_opts[0] = 10;
2389:     jac->as_amg_alpha_opts[1] = 1;
2390:     jac->as_amg_alpha_opts[2] = 6;
2391:     jac->as_amg_alpha_opts[3] = 6;
2392:     jac->as_amg_alpha_opts[4] = 4;
2393:     jac->as_amg_alpha_theta   = 0.25;
2394:     /* Scalar Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
2395:     jac->as_amg_beta_opts[0] = 10;
2396:     jac->as_amg_beta_opts[1] = 1;
2397:     jac->as_amg_beta_opts[2] = 6;
2398:     jac->as_amg_beta_opts[3] = 6;
2399:     jac->as_amg_beta_opts[4] = 4;
2400:     jac->as_amg_beta_theta   = 0.25;
2401:     PetscCallExternal(HYPRE_AMSSetPrintLevel, jac->hsolver, jac->as_print);
2402:     PetscCallExternal(HYPRE_AMSSetMaxIter, jac->hsolver, jac->as_max_iter);
2403:     PetscCallExternal(HYPRE_AMSSetCycleType, jac->hsolver, jac->ams_cycle_type);
2404:     PetscCallExternal(HYPRE_AMSSetTol, jac->hsolver, jac->as_tol);
2405:     PetscCallExternal(HYPRE_AMSSetSmoothingOptions, jac->hsolver, jac->as_relax_type, jac->as_relax_times, jac->as_relax_weight, jac->as_omega);
2406:     PetscCallExternal(HYPRE_AMSSetAlphaAMGOptions, jac->hsolver, jac->as_amg_alpha_opts[0], /* AMG coarsen type */
2407:                       jac->as_amg_alpha_opts[1],                                            /* AMG agg_levels */
2408:                       jac->as_amg_alpha_opts[2],                                            /* AMG relax_type */
2409:                       jac->as_amg_alpha_theta, jac->as_amg_alpha_opts[3],                   /* AMG interp_type */
2410:                       jac->as_amg_alpha_opts[4]);                                           /* AMG Pmax */
2411:     PetscCallExternal(HYPRE_AMSSetBetaAMGOptions, jac->hsolver, jac->as_amg_beta_opts[0],   /* AMG coarsen type */
2412:                       jac->as_amg_beta_opts[1],                                             /* AMG agg_levels */
2413:                       jac->as_amg_beta_opts[2],                                             /* AMG relax_type */
2414:                       jac->as_amg_beta_theta, jac->as_amg_beta_opts[3],                     /* AMG interp_type */
2415:                       jac->as_amg_beta_opts[4]);                                            /* AMG Pmax */
2416:     /* Zero conductivity */
2417:     jac->ams_beta_is_zero      = PETSC_FALSE;
2418:     jac->ams_beta_is_zero_part = PETSC_FALSE;
2419:     PetscFunctionReturn(PETSC_SUCCESS);
2420:   }
2421:   PetscCall(PetscStrcmp("ads", jac->hypre_type, &flag));
2422:   if (flag) {
2423:     PetscCallExternal(HYPRE_ADSCreate, &jac->hsolver);
2424:     pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ADS;
2425:     pc->ops->view           = PCView_HYPRE_ADS;
2426:     jac->destroy            = HYPRE_ADSDestroy;
2427:     jac->setup              = HYPRE_ADSSetup;
2428:     jac->solve              = HYPRE_ADSSolve;
2429:     jac->coords[0]          = NULL;
2430:     jac->coords[1]          = NULL;
2431:     jac->coords[2]          = NULL;
2432:     /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE ADS */
2433:     jac->as_print       = 0;
2434:     jac->as_max_iter    = 1;  /* used as a preconditioner */
2435:     jac->as_tol         = 0.; /* used as a preconditioner */
2436:     jac->ads_cycle_type = 13;
2437:     /* Smoothing options */
2438:     jac->as_relax_type   = 2;
2439:     jac->as_relax_times  = 1;
2440:     jac->as_relax_weight = 1.0;
2441:     jac->as_omega        = 1.0;
2442:     /* AMS solver parameters: cycle_type, coarsen type, agg_levels, relax_type, interp_type, Pmax */
2443:     jac->ams_cycle_type       = 14;
2444:     jac->as_amg_alpha_opts[0] = 10;
2445:     jac->as_amg_alpha_opts[1] = 1;
2446:     jac->as_amg_alpha_opts[2] = 6;
2447:     jac->as_amg_alpha_opts[3] = 6;
2448:     jac->as_amg_alpha_opts[4] = 4;
2449:     jac->as_amg_alpha_theta   = 0.25;
2450:     /* Vector Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */
2451:     jac->as_amg_beta_opts[0] = 10;
2452:     jac->as_amg_beta_opts[1] = 1;
2453:     jac->as_amg_beta_opts[2] = 6;
2454:     jac->as_amg_beta_opts[3] = 6;
2455:     jac->as_amg_beta_opts[4] = 4;
2456:     jac->as_amg_beta_theta   = 0.25;
2457:     PetscCallExternal(HYPRE_ADSSetPrintLevel, jac->hsolver, jac->as_print);
2458:     PetscCallExternal(HYPRE_ADSSetMaxIter, jac->hsolver, jac->as_max_iter);
2459:     PetscCallExternal(HYPRE_ADSSetCycleType, jac->hsolver, jac->ams_cycle_type);
2460:     PetscCallExternal(HYPRE_ADSSetTol, jac->hsolver, jac->as_tol);
2461:     PetscCallExternal(HYPRE_ADSSetSmoothingOptions, jac->hsolver, jac->as_relax_type, jac->as_relax_times, jac->as_relax_weight, jac->as_omega);
2462:     PetscCallExternal(HYPRE_ADSSetAMSOptions, jac->hsolver, jac->ams_cycle_type,      /* AMG coarsen type */
2463:                       jac->as_amg_alpha_opts[0],                                      /* AMG coarsen type */
2464:                       jac->as_amg_alpha_opts[1],                                      /* AMG agg_levels */
2465:                       jac->as_amg_alpha_opts[2],                                      /* AMG relax_type */
2466:                       jac->as_amg_alpha_theta, jac->as_amg_alpha_opts[3],             /* AMG interp_type */
2467:                       jac->as_amg_alpha_opts[4]);                                     /* AMG Pmax */
2468:     PetscCallExternal(HYPRE_ADSSetAMGOptions, jac->hsolver, jac->as_amg_beta_opts[0], /* AMG coarsen type */
2469:                       jac->as_amg_beta_opts[1],                                       /* AMG agg_levels */
2470:                       jac->as_amg_beta_opts[2],                                       /* AMG relax_type */
2471:                       jac->as_amg_beta_theta, jac->as_amg_beta_opts[3],               /* AMG interp_type */
2472:                       jac->as_amg_beta_opts[4]);                                      /* AMG Pmax */
2473:     PetscFunctionReturn(PETSC_SUCCESS);
2474:   }
2475:   PetscCall(PetscFree(jac->hypre_type));

2477:   jac->hypre_type = NULL;
2478:   SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown HYPRE preconditioner %s; Choices are euclid, ilu, pilut, parasails, boomeramg, ams, ads", name);
2479: }

2481: /*
2482:     It only gets here if the HYPRE type has not been set before the call to
2483:    ...SetFromOptions() which actually is most of the time
2484: */
2485: static PetscErrorCode PCSetFromOptions_HYPRE(PC pc, PetscOptionItems PetscOptionsObject)
2486: {
2487:   PetscInt    indx;
2488:   const char *type[] = {"ilu", "euclid", "pilut", "parasails", "boomeramg", "ams", "ads"};
2489:   PetscBool   flg;
2490:   PC_HYPRE   *jac = (PC_HYPRE *)pc->data;

2492:   PetscFunctionBegin;
2493:   PetscOptionsHeadBegin(PetscOptionsObject, "HYPRE preconditioner options");
2494:   PetscCall(PetscOptionsEList("-pc_hypre_type", "HYPRE preconditioner type", "PCHYPRESetType", type, PETSC_STATIC_ARRAY_LENGTH(type), "boomeramg", &indx, &flg));
2495:   if (flg) PetscCall(PCHYPRESetType_HYPRE(pc, type[indx]));
2496:   /*
2497:     Set the type if it was never set.
2498:   */
2499:   if (!jac->hypre_type) PetscCall(PCHYPRESetType_HYPRE(pc, "boomeramg"));
2500:   PetscTryTypeMethod(pc, setfromoptions, PetscOptionsObject);
2501:   PetscOptionsHeadEnd();
2502:   PetscFunctionReturn(PETSC_SUCCESS);
2503: }

2505: /*@
2506:   PCHYPRESetType - Sets which hypre preconditioner you wish to use

2508:   Input Parameters:
2509: + pc   - the preconditioner context
2510: - name - either euclid, ilu, pilut, parasails, boomeramg, ams, ads

2512:   Options Database Key:
2513: . pc_hypre_type - One of euclid, ilu, pilut, parasails, boomeramg, ams, ads

2515:   Level: intermediate

2517: .seealso: [](ch_ksp), `PCCreate()`, `PCSetType()`, `PCType`, `PC`, `PCHYPRE`
2518: @*/
2519: PetscErrorCode PCHYPRESetType(PC pc, const char name[])
2520: {
2521:   PetscFunctionBegin;
2523:   PetscAssertPointer(name, 2);
2524:   PetscTryMethod(pc, "PCHYPRESetType_C", (PC, const char[]), (pc, name));
2525:   PetscFunctionReturn(PETSC_SUCCESS);
2526: }

2528: /*@C
2529:   PCHYPREGetCFMarkers - Gets CF marker arrays for all levels (except the finest level)

2531:   Logically Collective

2533:   Input Parameter:
2534: . pc - the preconditioner context

2536:   Output Parameters:
2537: + n_per_level - the number of nodes per level (size of `num_levels`)
2538: - CFMarkers   - the Coarse/Fine Boolean arrays (size of `num_levels` - 1)

2540:   Note:
2541:   Caller is responsible for memory management of `n_per_level` and `CFMarkers` pointers. That is they should free them with `PetscFree()` when no longer needed.

2543:   Level: advanced

2545: .seealso: [](ch_ksp), `PC`, `PCMG`, `PCMGGetRestriction()`, `PCMGSetInterpolation()`, `PCMGGetRScale()`, `PCMGGetInterpolation()`, `PCGetInterpolations()`
2546: @*/
2547: PetscErrorCode PCHYPREGetCFMarkers(PC pc, PetscInt *n_per_level[], PetscBT *CFMarkers[])
2548: {
2549:   PetscFunctionBegin;
2551:   PetscAssertPointer(n_per_level, 2);
2552:   PetscAssertPointer(CFMarkers, 3);
2553:   PetscUseMethod(pc, "PCHYPREGetCFMarkers_C", (PC, PetscInt *[], PetscBT *[]), (pc, n_per_level, CFMarkers));
2554:   PetscFunctionReturn(PETSC_SUCCESS);
2555: }

2557: /*@
2558:   PCHYPREGetType - Gets which hypre preconditioner you are using

2560:   Input Parameter:
2561: . pc - the preconditioner context

2563:   Output Parameter:
2564: . name - either euclid, ilu, pilut, parasails, boomeramg, ams, ads

2566:   Level: intermediate

2568: .seealso: [](ch_ksp), `PCCreate()`, `PCHYPRESetType()`, `PCType`, `PC`, `PCHYPRE`
2569: @*/
2570: PetscErrorCode PCHYPREGetType(PC pc, const char *name[])
2571: {
2572:   PetscFunctionBegin;
2574:   PetscAssertPointer(name, 2);
2575:   PetscTryMethod(pc, "PCHYPREGetType_C", (PC, const char *[]), (pc, name));
2576:   PetscFunctionReturn(PETSC_SUCCESS);
2577: }

2579: /*@
2580:   PCMGGalerkinSetMatProductAlgorithm - Set type of SpGEMM for hypre to use on GPUs

2582:   Logically Collective

2584:   Input Parameters:
2585: + pc   - the hypre context
2586: - name - one of 'cusparse', 'hypre'

2588:   Options Database Key:
2589: . -pc_mg_galerkin_mat_product_algorithm <cusparse,hypre> - Type of SpGEMM to use in hypre

2591:   Level: intermediate

2593:   Developer Notes:
2594:   How the name starts with `PCMG`, should it not be `PCHYPREBoomerAMG`?

2596: .seealso: [](ch_ksp), `PCHYPRE`, `PCMGGalerkinGetMatProductAlgorithm()`
2597: @*/
2598: PetscErrorCode PCMGGalerkinSetMatProductAlgorithm(PC pc, const char name[])
2599: {
2600:   PetscFunctionBegin;
2602:   PetscTryMethod(pc, "PCMGGalerkinSetMatProductAlgorithm_C", (PC, const char[]), (pc, name));
2603:   PetscFunctionReturn(PETSC_SUCCESS);
2604: }

2606: /*@
2607:   PCMGGalerkinGetMatProductAlgorithm - Get type of SpGEMM for hypre to use on GPUs

2609:   Not Collective

2611:   Input Parameter:
2612: . pc - the multigrid context

2614:   Output Parameter:
2615: . name - one of 'cusparse', 'hypre'

2617:   Level: intermediate

2619: .seealso: [](ch_ksp), `PCHYPRE`, `PCMGGalerkinSetMatProductAlgorithm()`
2620: @*/
2621: PetscErrorCode PCMGGalerkinGetMatProductAlgorithm(PC pc, const char *name[])
2622: {
2623:   PetscFunctionBegin;
2625:   PetscTryMethod(pc, "PCMGGalerkinGetMatProductAlgorithm_C", (PC, const char *[]), (pc, name));
2626:   PetscFunctionReturn(PETSC_SUCCESS);
2627: }

2629: /*MC
2630:      PCHYPRE - Allows you to use the matrix element based preconditioners in the LLNL package hypre as PETSc `PC`

2632:    Options Database Keys:
2633: +   -pc_hypre_type - One of `euclid`, `pilut`, `parasails`, `boomeramg`, `ams`, or `ads`
2634: .   -pc_hypre_boomeramg_nodal_coarsen <n> - where n is from 1 to 6 (see `HYPRE_BoomerAMGSetNodal()`)
2635: .   -pc_hypre_boomeramg_vec_interp_variant <v> - where v is from 1 to 3 (see `HYPRE_BoomerAMGSetInterpVecVariant()`)
2636: -   Many others - run with `-pc_type hypre` `-pc_hypre_type XXX` `-help` to see options for the XXX preconditioner

2638:    Level: intermediate

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

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

2653:            Note that the option `-pc_hypre_boomeramg_relax_type_all` defaults to symmetric relaxation
2654:            (symmetric-SOR/Jacobi), which is required for Krylov solvers like CG that expect symmetry.
2655:            Otherwise, you may want to use `-pc_hypre_boomeramg_relax_type_all SOR/Jacobi`.

2657:           `MatSetNearNullSpace()` - if you provide a near null space to your matrix it is ignored by hypre UNLESS you also use
2658:           the following two options: `-pc_hypre_boomeramg_nodal_coarsen <n> -pc_hypre_boomeramg_vec_interp_variant <v>`

2660:           See `PCPFMG`, `PCSMG`, and `PCSYSPFMG` for access to hypre's other (nonalgebraic) multigrid solvers

2662:           For `PCHYPRE` type of `ams` or `ads` auxiliary data must be provided to the preconditioner with `PCHYPRESetDiscreteGradient()`,
2663:           `PCHYPRESetDiscreteCurl()`, `PCHYPRESetInterpolations()`, `PCHYPRESetAlphaPoissonMatrix()`, `PCHYPRESetBetaPoissonMatrix()`, `PCHYPRESetEdgeConstantVectors()`,
2664:           `PCHYPREAMSSetInteriorNodes()`

2666:   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
2667:   since it is usually slower, one should use a `KSPType` of `KSPRICHARDSON`
2668:   (or equivalently `-ksp_type richardson`) to achieve this. Using `KSPPREONLY` will not work since it only applies a single cycle of multigrid.

2670:    PETSc provides its own geometric and algebraic multigrid solvers `PCMG` and `PCGAMG`, also see `PCHMG` which is useful for certain multicomponent problems

2672:    GPU Notes:
2673:      To configure hypre BoomerAMG so that it can utilize NVIDIA GPUs run ./configure --download-hypre --with-cuda
2674:      Then pass `VECCUDA` vectors and `MATAIJCUSPARSE` matrices to the solvers and PETSc will automatically utilize hypre's GPU solvers.

2676:      To configure hypre BoomerAMG so that it can utilize AMD GPUs run ./configure --download-hypre --with-hip
2677:      Then pass `VECHIP` vectors to the solvers and PETSc will automatically utilize hypre's GPU solvers.

2679: .seealso: [](ch_ksp), `PCCreate()`, `PCSetType()`, `PCType`, `PC`, `PCHYPRESetType()`, `PCPFMG`, `PCGAMG`, `PCSYSPFMG`, `PCSMG`, `PCHYPRESetDiscreteGradient()`,
2680:           `PCHYPRESetDiscreteCurl()`, `PCHYPRESetInterpolations()`, `PCHYPRESetAlphaPoissonMatrix()`, `PCHYPRESetBetaPoissonMatrix()`, `PCHYPRESetEdgeConstantVectors()`,
2681:           PCHYPREAMSSetInteriorNodes()
2682: M*/

2684: PETSC_EXTERN PetscErrorCode PCCreate_HYPRE(PC pc)
2685: {
2686:   PC_HYPRE *jac;

2688:   PetscFunctionBegin;
2689:   PetscCall(PetscNew(&jac));

2691:   pc->data                = jac;
2692:   pc->ops->reset          = PCReset_HYPRE;
2693:   pc->ops->destroy        = PCDestroy_HYPRE;
2694:   pc->ops->setfromoptions = PCSetFromOptions_HYPRE;
2695:   pc->ops->setup          = PCSetUp_HYPRE;
2696:   pc->ops->apply          = PCApply_HYPRE;
2697:   jac->hypre_type         = NULL;
2698:   jac->comm_hypre         = MPI_COMM_NULL;
2699:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetType_C", PCHYPRESetType_HYPRE));
2700:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPREGetType_C", PCHYPREGetType_HYPRE));
2701:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSetCoordinates_C", PCSetCoordinates_HYPRE));
2702:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetDiscreteGradient_C", PCHYPRESetDiscreteGradient_HYPRE));
2703:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetDiscreteCurl_C", PCHYPRESetDiscreteCurl_HYPRE));
2704:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetInterpolations_C", PCHYPRESetInterpolations_HYPRE));
2705:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetEdgeConstantVectors_C", PCHYPRESetEdgeConstantVectors_HYPRE));
2706:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPREAMSSetInteriorNodes_C", PCHYPREAMSSetInteriorNodes_HYPRE));
2707:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHYPRESetPoissonMatrix_C", PCHYPRESetPoissonMatrix_HYPRE));
2708:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCMGGalerkinSetMatProductAlgorithm_C", PCMGGalerkinSetMatProductAlgorithm_HYPRE_BoomerAMG));
2709:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCMGGalerkinGetMatProductAlgorithm_C", PCMGGalerkinGetMatProductAlgorithm_HYPRE_BoomerAMG));
2710: #if defined(PETSC_HAVE_HYPRE_DEVICE)
2711:   #if defined(HYPRE_USING_HIP)
2712:   PetscCall(PetscDeviceInitialize(PETSC_DEVICE_HIP));
2713:   #endif
2714:   #if defined(HYPRE_USING_CUDA)
2715:   PetscCall(PetscDeviceInitialize(PETSC_DEVICE_CUDA));
2716:   #endif
2717: #endif
2718:   PetscHYPREInitialize();
2719:   PetscFunctionReturn(PETSC_SUCCESS);
2720: }

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

2726:   /* keep copy of PFMG options used so may view them */
2727:   PetscInt  its;
2728:   PetscReal tol;
2729:   PetscInt  relax_type;
2730:   PetscInt  rap_type;
2731:   PetscInt  num_pre_relax, num_post_relax;
2732:   PetscInt  max_levels;
2733:   PetscInt  skip_relax;
2734:   PetscBool print_statistics;
2735: } PC_PFMG;

2737: static PetscErrorCode PCDestroy_PFMG(PC pc)
2738: {
2739:   PC_PFMG *ex = (PC_PFMG *)pc->data;

2741:   PetscFunctionBegin;
2742:   if (ex->hsolver) PetscCallExternal(HYPRE_StructPFMGDestroy, ex->hsolver);
2743:   PetscCall(PetscCommRestoreComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
2744:   PetscCall(PetscFree(pc->data));
2745:   PetscFunctionReturn(PETSC_SUCCESS);
2746: }

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

2751: static PetscErrorCode PCView_PFMG(PC pc, PetscViewer viewer)
2752: {
2753:   PetscBool isascii;
2754:   PC_PFMG  *ex = (PC_PFMG *)pc->data;

2756:   PetscFunctionBegin;
2757:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
2758:   if (isascii) {
2759:     PetscCall(PetscViewerASCIIPrintf(viewer, "  HYPRE PFMG preconditioning\n"));
2760:     PetscCall(PetscViewerASCIIPrintf(viewer, "    max iterations %" PetscInt_FMT "\n", ex->its));
2761:     PetscCall(PetscViewerASCIIPrintf(viewer, "    tolerance %g\n", ex->tol));
2762:     PetscCall(PetscViewerASCIIPrintf(viewer, "    relax type %s\n", PFMGRelaxType[ex->relax_type]));
2763:     PetscCall(PetscViewerASCIIPrintf(viewer, "    RAP type %s\n", PFMGRAPType[ex->rap_type]));
2764:     PetscCall(PetscViewerASCIIPrintf(viewer, "    number pre-relax %" PetscInt_FMT " post-relax %" PetscInt_FMT "\n", ex->num_pre_relax, ex->num_post_relax));
2765:     PetscCall(PetscViewerASCIIPrintf(viewer, "    max levels %" PetscInt_FMT "\n", ex->max_levels));
2766:     PetscCall(PetscViewerASCIIPrintf(viewer, "    skip relax %" PetscInt_FMT "\n", ex->skip_relax));
2767:   }
2768:   PetscFunctionReturn(PETSC_SUCCESS);
2769: }

2771: static PetscErrorCode PCSetFromOptions_PFMG(PC pc, PetscOptionItems PetscOptionsObject)
2772: {
2773:   PC_PFMG *ex = (PC_PFMG *)pc->data;

2775:   PetscFunctionBegin;
2776:   PetscOptionsHeadBegin(PetscOptionsObject, "PFMG options");
2777:   PetscCall(PetscOptionsBool("-pc_pfmg_print_statistics", "Print statistics", "HYPRE_StructPFMGSetPrintLevel", ex->print_statistics, &ex->print_statistics, NULL));
2778:   PetscCall(PetscOptionsInt("-pc_pfmg_its", "Number of iterations of PFMG to use as preconditioner", "HYPRE_StructPFMGSetMaxIter", ex->its, &ex->its, NULL));
2779:   PetscCallExternal(HYPRE_StructPFMGSetMaxIter, ex->hsolver, ex->its);
2780:   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));
2781:   PetscCallExternal(HYPRE_StructPFMGSetNumPreRelax, ex->hsolver, ex->num_pre_relax);
2782:   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));
2783:   PetscCallExternal(HYPRE_StructPFMGSetNumPostRelax, ex->hsolver, ex->num_post_relax);

2785:   PetscCall(PetscOptionsInt("-pc_pfmg_max_levels", "Max Levels for MG hierarchy", "HYPRE_StructPFMGSetMaxLevels", ex->max_levels, &ex->max_levels, NULL));
2786:   PetscCallExternal(HYPRE_StructPFMGSetMaxLevels, ex->hsolver, ex->max_levels);

2788:   PetscCall(PetscOptionsReal("-pc_pfmg_tol", "Tolerance of PFMG", "HYPRE_StructPFMGSetTol", ex->tol, &ex->tol, NULL));
2789:   PetscCallExternal(HYPRE_StructPFMGSetTol, ex->hsolver, ex->tol);
2790:   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));
2791:   PetscCallExternal(HYPRE_StructPFMGSetRelaxType, ex->hsolver, ex->relax_type);
2792:   PetscCall(PetscOptionsEList("-pc_pfmg_rap_type", "RAP type", "HYPRE_StructPFMGSetRAPType", PFMGRAPType, PETSC_STATIC_ARRAY_LENGTH(PFMGRAPType), PFMGRAPType[ex->rap_type], &ex->rap_type, NULL));
2793:   PetscCallExternal(HYPRE_StructPFMGSetRAPType, ex->hsolver, ex->rap_type);
2794:   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));
2795:   PetscCallExternal(HYPRE_StructPFMGSetSkipRelax, ex->hsolver, ex->skip_relax);
2796:   PetscOptionsHeadEnd();
2797:   PetscFunctionReturn(PETSC_SUCCESS);
2798: }

2800: static PetscErrorCode PCApply_PFMG(PC pc, Vec x, Vec y)
2801: {
2802:   PC_PFMG           *ex = (PC_PFMG *)pc->data;
2803:   PetscScalar       *yy;
2804:   const PetscScalar *xx;
2805:   PetscInt           ilower[3], iupper[3];
2806:   HYPRE_Int          hlower[3], hupper[3];
2807:   Mat_HYPREStruct   *mx = (Mat_HYPREStruct *)pc->pmat->data;

2809:   PetscFunctionBegin;
2810:   PetscCall(PetscCitationsRegister(hypreCitation, &cite));
2811:   PetscCall(DMDAGetCorners(mx->da, &ilower[0], &ilower[1], &ilower[2], &iupper[0], &iupper[1], &iupper[2]));
2812:   /* when HYPRE_MIXEDINT is defined, sizeof(HYPRE_Int) == 32 */
2813:   iupper[0] += ilower[0] - 1;
2814:   iupper[1] += ilower[1] - 1;
2815:   iupper[2] += ilower[2] - 1;
2816:   hlower[0] = (HYPRE_Int)ilower[0];
2817:   hlower[1] = (HYPRE_Int)ilower[1];
2818:   hlower[2] = (HYPRE_Int)ilower[2];
2819:   hupper[0] = (HYPRE_Int)iupper[0];
2820:   hupper[1] = (HYPRE_Int)iupper[1];
2821:   hupper[2] = (HYPRE_Int)iupper[2];

2823:   /* copy x values over to hypre */
2824:   PetscCallExternal(HYPRE_StructVectorSetConstantValues, mx->hb, 0.0);
2825:   PetscCall(VecGetArrayRead(x, &xx));
2826:   PetscCallExternal(HYPRE_StructVectorSetBoxValues, mx->hb, hlower, hupper, (HYPRE_Complex *)xx);
2827:   PetscCall(VecRestoreArrayRead(x, &xx));
2828:   PetscCallExternal(HYPRE_StructVectorAssemble, mx->hb);
2829:   PetscCallExternal(HYPRE_StructPFMGSolve, ex->hsolver, mx->hmat, mx->hb, mx->hx);

2831:   /* copy solution values back to PETSc */
2832:   PetscCall(VecGetArray(y, &yy));
2833:   PetscCallExternal(HYPRE_StructVectorGetBoxValues, mx->hx, hlower, hupper, (HYPRE_Complex *)yy);
2834:   PetscCall(VecRestoreArray(y, &yy));
2835:   PetscFunctionReturn(PETSC_SUCCESS);
2836: }

2838: 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)
2839: {
2840:   PC_PFMG  *jac = (PC_PFMG *)pc->data;
2841:   HYPRE_Int oits;

2843:   PetscFunctionBegin;
2844:   PetscCall(PetscCitationsRegister(hypreCitation, &cite));
2845:   PetscCallExternal(HYPRE_StructPFMGSetMaxIter, jac->hsolver, its * jac->its);
2846:   PetscCallExternal(HYPRE_StructPFMGSetTol, jac->hsolver, rtol);

2848:   PetscCall(PCApply_PFMG(pc, b, y));
2849:   PetscCallExternal(HYPRE_StructPFMGGetNumIterations, jac->hsolver, &oits);
2850:   *outits = oits;
2851:   if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
2852:   else *reason = PCRICHARDSON_CONVERGED_RTOL;
2853:   PetscCallExternal(HYPRE_StructPFMGSetTol, jac->hsolver, jac->tol);
2854:   PetscCallExternal(HYPRE_StructPFMGSetMaxIter, jac->hsolver, jac->its);
2855:   PetscFunctionReturn(PETSC_SUCCESS);
2856: }

2858: static PetscErrorCode PCSetUp_PFMG(PC pc)
2859: {
2860:   PC_PFMG         *ex = (PC_PFMG *)pc->data;
2861:   Mat_HYPREStruct *mx = (Mat_HYPREStruct *)pc->pmat->data;
2862:   PetscBool        flg;

2864:   PetscFunctionBegin;
2865:   PetscCall(PetscObjectTypeCompare((PetscObject)pc->pmat, MATHYPRESTRUCT, &flg));
2866:   PetscCheck(flg, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_INCOMP, "Must use MATHYPRESTRUCT with this preconditioner");

2868:   /* create the hypre solver object and set its information */
2869:   if (ex->hsolver) PetscCallExternal(HYPRE_StructPFMGDestroy, ex->hsolver);
2870:   PetscCallExternal(HYPRE_StructPFMGCreate, ex->hcomm, &ex->hsolver);

2872:   // Print Hypre statistics about the solve process
2873:   if (ex->print_statistics) PetscCallExternal(HYPRE_StructPFMGSetPrintLevel, ex->hsolver, 3);

2875:   // The hypre options must be repeated here because the StructPFMG was destroyed and recreated
2876:   PetscCallExternal(HYPRE_StructPFMGSetMaxIter, ex->hsolver, ex->its);
2877:   PetscCallExternal(HYPRE_StructPFMGSetNumPreRelax, ex->hsolver, ex->num_pre_relax);
2878:   PetscCallExternal(HYPRE_StructPFMGSetNumPostRelax, ex->hsolver, ex->num_post_relax);
2879:   PetscCallExternal(HYPRE_StructPFMGSetMaxLevels, ex->hsolver, ex->max_levels);
2880:   PetscCallExternal(HYPRE_StructPFMGSetTol, ex->hsolver, ex->tol);
2881:   PetscCallExternal(HYPRE_StructPFMGSetRelaxType, ex->hsolver, ex->relax_type);
2882:   PetscCallExternal(HYPRE_StructPFMGSetRAPType, ex->hsolver, ex->rap_type);

2884:   PetscCallExternal(HYPRE_StructPFMGSetup, ex->hsolver, mx->hmat, mx->hb, mx->hx);
2885:   PetscCallExternal(HYPRE_StructPFMGSetZeroGuess, ex->hsolver);
2886:   PetscFunctionReturn(PETSC_SUCCESS);
2887: }

2889: /*MC
2890:      PCPFMG - the hypre PFMG multigrid solver

2892:    Options Database Keys:
2893: + -pc_pfmg_its <its> - number of iterations of PFMG to use as preconditioner
2894: . -pc_pfmg_num_pre_relax <steps> - number of smoothing steps before coarse grid solve
2895: . -pc_pfmg_num_post_relax <steps> - number of smoothing steps after coarse grid solve
2896: . -pc_pfmg_tol <tol> - tolerance of PFMG
2897: . -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
2898: . -pc_pfmg_rap_type - type of coarse matrix generation, one of Galerkin,non-Galerkin
2899: - -pc_pfmg_skip_relax - skip relaxation on certain grids for isotropic problems. This can greatly improve efficiency by eliminating unnecessary relaxations
2900:                         when the underlying problem is isotropic, one of 0,1

2902:    Level: advanced

2904:    Notes:
2905:    This is for CELL-centered descretizations

2907:    See `PCSYSPFMG` for a version suitable for systems of PDEs, and `PCSMG`

2909:    See `PCHYPRE` for hypre's BoomerAMG algebraic multigrid solver

2911:    This must be used with the `MATHYPRESTRUCT` matrix type.

2913:    This provides only some of the functionality of PFMG, it supports only one block per process defined by a PETSc `DMDA`.

2915: .seealso: [](ch_ksp), `PCMG`, `MATHYPRESTRUCT`, `PCHYPRE`, `PCGAMG`, `PCSYSPFMG`, `PCSMG`
2916: M*/

2918: PETSC_EXTERN PetscErrorCode PCCreate_PFMG(PC pc)
2919: {
2920:   PC_PFMG *ex;

2922:   PetscFunctionBegin;
2923:   PetscCall(PetscNew(&ex));
2924:   pc->data = ex;

2926:   ex->its              = 1;
2927:   ex->tol              = 1.e-8;
2928:   ex->relax_type       = 1;
2929:   ex->rap_type         = 0;
2930:   ex->num_pre_relax    = 1;
2931:   ex->num_post_relax   = 1;
2932:   ex->max_levels       = 0;
2933:   ex->skip_relax       = 0;
2934:   ex->print_statistics = PETSC_FALSE;

2936:   pc->ops->setfromoptions  = PCSetFromOptions_PFMG;
2937:   pc->ops->view            = PCView_PFMG;
2938:   pc->ops->destroy         = PCDestroy_PFMG;
2939:   pc->ops->apply           = PCApply_PFMG;
2940:   pc->ops->applyrichardson = PCApplyRichardson_PFMG;
2941:   pc->ops->setup           = PCSetUp_PFMG;

2943:   PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
2944:   PetscHYPREInitialize();
2945:   PetscCallExternal(HYPRE_StructPFMGCreate, ex->hcomm, &ex->hsolver);
2946:   PetscFunctionReturn(PETSC_SUCCESS);
2947: }

2949: /* we know we are working with a HYPRE_SStructMatrix */
2950: typedef struct {
2951:   MPI_Comm            hcomm; /* does not share comm with HYPRE_SStructMatrix because need to create solver before getting matrix */
2952:   HYPRE_SStructSolver ss_solver;

2954:   /* keep copy of SYSPFMG options used so may view them */
2955:   PetscInt  its;
2956:   PetscReal tol;
2957:   PetscInt  relax_type;
2958:   PetscInt  num_pre_relax, num_post_relax;
2959: } PC_SysPFMG;

2961: static PetscErrorCode PCDestroy_SysPFMG(PC pc)
2962: {
2963:   PC_SysPFMG *ex = (PC_SysPFMG *)pc->data;

2965:   PetscFunctionBegin;
2966:   if (ex->ss_solver) PetscCallExternal(HYPRE_SStructSysPFMGDestroy, ex->ss_solver);
2967:   PetscCall(PetscCommRestoreComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
2968:   PetscCall(PetscFree(pc->data));
2969:   PetscFunctionReturn(PETSC_SUCCESS);
2970: }

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

2974: static PetscErrorCode PCView_SysPFMG(PC pc, PetscViewer viewer)
2975: {
2976:   PetscBool   isascii;
2977:   PC_SysPFMG *ex = (PC_SysPFMG *)pc->data;

2979:   PetscFunctionBegin;
2980:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
2981:   if (isascii) {
2982:     PetscCall(PetscViewerASCIIPrintf(viewer, "  HYPRE SysPFMG preconditioning\n"));
2983:     PetscCall(PetscViewerASCIIPrintf(viewer, "  max iterations %" PetscInt_FMT "\n", ex->its));
2984:     PetscCall(PetscViewerASCIIPrintf(viewer, "  tolerance %g\n", ex->tol));
2985:     PetscCall(PetscViewerASCIIPrintf(viewer, "  relax type %s\n", PFMGRelaxType[ex->relax_type]));
2986:     PetscCall(PetscViewerASCIIPrintf(viewer, "  number pre-relax %" PetscInt_FMT " post-relax %" PetscInt_FMT "\n", ex->num_pre_relax, ex->num_post_relax));
2987:   }
2988:   PetscFunctionReturn(PETSC_SUCCESS);
2989: }

2991: static PetscErrorCode PCSetFromOptions_SysPFMG(PC pc, PetscOptionItems PetscOptionsObject)
2992: {
2993:   PC_SysPFMG *ex  = (PC_SysPFMG *)pc->data;
2994:   PetscBool   flg = PETSC_FALSE;

2996:   PetscFunctionBegin;
2997:   PetscOptionsHeadBegin(PetscOptionsObject, "SysPFMG options");
2998:   PetscCall(PetscOptionsBool("-pc_syspfmg_print_statistics", "Print statistics", "HYPRE_SStructSysPFMGSetPrintLevel", flg, &flg, NULL));
2999:   if (flg) PetscCallExternal(HYPRE_SStructSysPFMGSetPrintLevel, ex->ss_solver, 3);
3000:   PetscCall(PetscOptionsInt("-pc_syspfmg_its", "Number of iterations of SysPFMG to use as preconditioner", "HYPRE_SStructSysPFMGSetMaxIter", ex->its, &ex->its, NULL));
3001:   PetscCallExternal(HYPRE_SStructSysPFMGSetMaxIter, ex->ss_solver, ex->its);
3002:   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));
3003:   PetscCallExternal(HYPRE_SStructSysPFMGSetNumPreRelax, ex->ss_solver, ex->num_pre_relax);
3004:   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));
3005:   PetscCallExternal(HYPRE_SStructSysPFMGSetNumPostRelax, ex->ss_solver, ex->num_post_relax);

3007:   PetscCall(PetscOptionsReal("-pc_syspfmg_tol", "Tolerance of SysPFMG", "HYPRE_SStructSysPFMGSetTol", ex->tol, &ex->tol, NULL));
3008:   PetscCallExternal(HYPRE_SStructSysPFMGSetTol, ex->ss_solver, ex->tol);
3009:   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));
3010:   PetscCallExternal(HYPRE_SStructSysPFMGSetRelaxType, ex->ss_solver, ex->relax_type);
3011:   PetscOptionsHeadEnd();
3012:   PetscFunctionReturn(PETSC_SUCCESS);
3013: }

3015: static PetscErrorCode PCApply_SysPFMG(PC pc, Vec x, Vec y)
3016: {
3017:   PC_SysPFMG        *ex = (PC_SysPFMG *)pc->data;
3018:   PetscScalar       *yy;
3019:   const PetscScalar *xx;
3020:   PetscInt           ilower[3], iupper[3];
3021:   HYPRE_Int          hlower[3], hupper[3];
3022:   Mat_HYPRESStruct  *mx       = (Mat_HYPRESStruct *)pc->pmat->data;
3023:   PetscInt           ordering = mx->dofs_order;
3024:   PetscInt           nvars    = mx->nvars;
3025:   PetscInt           part     = 0;
3026:   PetscInt           size;
3027:   PetscInt           i;

3029:   PetscFunctionBegin;
3030:   PetscCall(PetscCitationsRegister(hypreCitation, &cite));
3031:   PetscCall(DMDAGetCorners(mx->da, &ilower[0], &ilower[1], &ilower[2], &iupper[0], &iupper[1], &iupper[2]));
3032:   /* when HYPRE_MIXEDINT is defined, sizeof(HYPRE_Int) == 32 */
3033:   iupper[0] += ilower[0] - 1;
3034:   iupper[1] += ilower[1] - 1;
3035:   iupper[2] += ilower[2] - 1;
3036:   hlower[0] = (HYPRE_Int)ilower[0];
3037:   hlower[1] = (HYPRE_Int)ilower[1];
3038:   hlower[2] = (HYPRE_Int)ilower[2];
3039:   hupper[0] = (HYPRE_Int)iupper[0];
3040:   hupper[1] = (HYPRE_Int)iupper[1];
3041:   hupper[2] = (HYPRE_Int)iupper[2];

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

3046:   /* copy x values over to hypre for variable ordering */
3047:   if (ordering) {
3048:     PetscCallExternal(HYPRE_SStructVectorSetConstantValues, mx->ss_b, 0.0);
3049:     PetscCall(VecGetArrayRead(x, &xx));
3050:     for (i = 0; i < nvars; i++) PetscCallExternal(HYPRE_SStructVectorSetBoxValues, mx->ss_b, part, hlower, hupper, i, (HYPRE_Complex *)(xx + (size * i)));
3051:     PetscCall(VecRestoreArrayRead(x, &xx));
3052:     PetscCallExternal(HYPRE_SStructVectorAssemble, mx->ss_b);
3053:     PetscCallExternal(HYPRE_SStructMatrixMatvec, 1.0, mx->ss_mat, mx->ss_b, 0.0, mx->ss_x);
3054:     PetscCallExternal(HYPRE_SStructSysPFMGSolve, ex->ss_solver, mx->ss_mat, mx->ss_b, mx->ss_x);

3056:     /* copy solution values back to PETSc */
3057:     PetscCall(VecGetArray(y, &yy));
3058:     for (i = 0; i < nvars; i++) PetscCallExternal(HYPRE_SStructVectorGetBoxValues, mx->ss_x, part, hlower, hupper, i, (HYPRE_Complex *)(yy + (size * i)));
3059:     PetscCall(VecRestoreArray(y, &yy));
3060:   } else { /* nodal ordering must be mapped to variable ordering for sys_pfmg */
3061:     PetscScalar *z;
3062:     PetscInt     j, k;

3064:     PetscCall(PetscMalloc1(nvars * size, &z));
3065:     PetscCallExternal(HYPRE_SStructVectorSetConstantValues, mx->ss_b, 0.0);
3066:     PetscCall(VecGetArrayRead(x, &xx));

3068:     /* transform nodal to hypre's variable ordering for sys_pfmg */
3069:     for (i = 0; i < size; i++) {
3070:       k = i * nvars;
3071:       for (j = 0; j < nvars; j++) z[j * size + i] = xx[k + j];
3072:     }
3073:     for (i = 0; i < nvars; i++) PetscCallExternal(HYPRE_SStructVectorSetBoxValues, mx->ss_b, part, hlower, hupper, i, (HYPRE_Complex *)(z + (size * i)));
3074:     PetscCall(VecRestoreArrayRead(x, &xx));
3075:     PetscCallExternal(HYPRE_SStructVectorAssemble, mx->ss_b);
3076:     PetscCallExternal(HYPRE_SStructSysPFMGSolve, ex->ss_solver, mx->ss_mat, mx->ss_b, mx->ss_x);

3078:     /* copy solution values back to PETSc */
3079:     PetscCall(VecGetArray(y, &yy));
3080:     for (i = 0; i < nvars; i++) PetscCallExternal(HYPRE_SStructVectorGetBoxValues, mx->ss_x, part, hlower, hupper, i, (HYPRE_Complex *)(z + (size * i)));
3081:     /* transform hypre's variable ordering for sys_pfmg to nodal ordering */
3082:     for (i = 0; i < size; i++) {
3083:       k = i * nvars;
3084:       for (j = 0; j < nvars; j++) yy[k + j] = z[j * size + i];
3085:     }
3086:     PetscCall(VecRestoreArray(y, &yy));
3087:     PetscCall(PetscFree(z));
3088:   }
3089:   PetscFunctionReturn(PETSC_SUCCESS);
3090: }

3092: 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)
3093: {
3094:   PC_SysPFMG *jac = (PC_SysPFMG *)pc->data;
3095:   HYPRE_Int   oits;

3097:   PetscFunctionBegin;
3098:   PetscCall(PetscCitationsRegister(hypreCitation, &cite));
3099:   PetscCallExternal(HYPRE_SStructSysPFMGSetMaxIter, jac->ss_solver, its * jac->its);
3100:   PetscCallExternal(HYPRE_SStructSysPFMGSetTol, jac->ss_solver, rtol);
3101:   PetscCall(PCApply_SysPFMG(pc, b, y));
3102:   PetscCallExternal(HYPRE_SStructSysPFMGGetNumIterations, jac->ss_solver, &oits);
3103:   *outits = oits;
3104:   if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
3105:   else *reason = PCRICHARDSON_CONVERGED_RTOL;
3106:   PetscCallExternal(HYPRE_SStructSysPFMGSetTol, jac->ss_solver, jac->tol);
3107:   PetscCallExternal(HYPRE_SStructSysPFMGSetMaxIter, jac->ss_solver, jac->its);
3108:   PetscFunctionReturn(PETSC_SUCCESS);
3109: }

3111: static PetscErrorCode PCSetUp_SysPFMG(PC pc)
3112: {
3113:   PC_SysPFMG       *ex = (PC_SysPFMG *)pc->data;
3114:   Mat_HYPRESStruct *mx = (Mat_HYPRESStruct *)pc->pmat->data;
3115:   PetscBool         flg;

3117:   PetscFunctionBegin;
3118:   PetscCall(PetscObjectTypeCompare((PetscObject)pc->pmat, MATHYPRESSTRUCT, &flg));
3119:   PetscCheck(flg, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_INCOMP, "Must use MATHYPRESSTRUCT with this preconditioner");

3121:   /* create the hypre sstruct solver object and set its information */
3122:   if (ex->ss_solver) PetscCallExternal(HYPRE_SStructSysPFMGDestroy, ex->ss_solver);
3123:   PetscCallExternal(HYPRE_SStructSysPFMGCreate, ex->hcomm, &ex->ss_solver);
3124:   PetscCallExternal(HYPRE_SStructSysPFMGSetZeroGuess, ex->ss_solver);
3125:   PetscCallExternal(HYPRE_SStructSysPFMGSetup, ex->ss_solver, mx->ss_mat, mx->ss_b, mx->ss_x);
3126:   PetscFunctionReturn(PETSC_SUCCESS);
3127: }

3129: /*MC
3130:      PCSYSPFMG - the hypre SysPFMG multigrid solver

3132:    Level: advanced

3134:    Options Database Keys:
3135: + -pc_syspfmg_its <its> - number of iterations of SysPFMG to use as preconditioner
3136: . -pc_syspfmg_num_pre_relax <steps> - number of smoothing steps before coarse grid
3137: . -pc_syspfmg_num_post_relax <steps> - number of smoothing steps after coarse grid
3138: . -pc_syspfmg_tol <tol> - tolerance of SysPFMG
3139: - -pc_syspfmg_relax_type <Weighted-Jacobi,Red/Black-Gauss-Seidel> - relaxation type for the up and down cycles

3141:    Notes:
3142:    See `PCPFMG` for hypre's PFMG that works for a scalar PDE and `PCSMG`

3144:    See `PCHYPRE` for hypre's BoomerAMG algebraic multigrid solver

3146:    This is for CELL-centered descretizations

3148:    This must be used with the `MATHYPRESSTRUCT` matrix type.

3150:    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`.

3152: .seealso: [](ch_ksp), `PCMG`, `MATHYPRESSTRUCT`, `PCPFMG`, `PCHYPRE`, `PCGAMG`, `PCSMG`
3153: M*/

3155: PETSC_EXTERN PetscErrorCode PCCreate_SysPFMG(PC pc)
3156: {
3157:   PC_SysPFMG *ex;

3159:   PetscFunctionBegin;
3160:   PetscCall(PetscNew(&ex));
3161:   pc->data = ex;

3163:   ex->its            = 1;
3164:   ex->tol            = 1.e-8;
3165:   ex->relax_type     = 1;
3166:   ex->num_pre_relax  = 1;
3167:   ex->num_post_relax = 1;

3169:   pc->ops->setfromoptions  = PCSetFromOptions_SysPFMG;
3170:   pc->ops->view            = PCView_SysPFMG;
3171:   pc->ops->destroy         = PCDestroy_SysPFMG;
3172:   pc->ops->apply           = PCApply_SysPFMG;
3173:   pc->ops->applyrichardson = PCApplyRichardson_SysPFMG;
3174:   pc->ops->setup           = PCSetUp_SysPFMG;

3176:   PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
3177:   PetscHYPREInitialize();
3178:   PetscCallExternal(HYPRE_SStructSysPFMGCreate, ex->hcomm, &ex->ss_solver);
3179:   PetscFunctionReturn(PETSC_SUCCESS);
3180: }

3182: /* PC SMG */
3183: typedef struct {
3184:   MPI_Comm           hcomm; /* does not share comm with HYPRE_StructMatrix because need to create solver before getting matrix */
3185:   HYPRE_StructSolver hsolver;
3186:   PetscInt           its; /* keep copy of SMG options used so may view them */
3187:   PetscReal          tol;
3188:   PetscBool          print_statistics;
3189:   PetscInt           num_pre_relax, num_post_relax;
3190: } PC_SMG;

3192: static PetscErrorCode PCDestroy_SMG(PC pc)
3193: {
3194:   PC_SMG *ex = (PC_SMG *)pc->data;

3196:   PetscFunctionBegin;
3197:   if (ex->hsolver) PetscCallExternal(HYPRE_StructSMGDestroy, ex->hsolver);
3198:   PetscCall(PetscCommRestoreComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
3199:   PetscCall(PetscFree(pc->data));
3200:   PetscFunctionReturn(PETSC_SUCCESS);
3201: }

3203: static PetscErrorCode PCView_SMG(PC pc, PetscViewer viewer)
3204: {
3205:   PetscBool isascii;
3206:   PC_SMG   *ex = (PC_SMG *)pc->data;

3208:   PetscFunctionBegin;
3209:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
3210:   if (isascii) {
3211:     PetscCall(PetscViewerASCIIPrintf(viewer, "  HYPRE SMG preconditioning\n"));
3212:     PetscCall(PetscViewerASCIIPrintf(viewer, "    max iterations %" PetscInt_FMT "\n", ex->its));
3213:     PetscCall(PetscViewerASCIIPrintf(viewer, "    tolerance %g\n", ex->tol));
3214:     PetscCall(PetscViewerASCIIPrintf(viewer, "    number pre-relax %" PetscInt_FMT " post-relax %" PetscInt_FMT "\n", ex->num_pre_relax, ex->num_post_relax));
3215:   }
3216:   PetscFunctionReturn(PETSC_SUCCESS);
3217: }

3219: static PetscErrorCode PCSetFromOptions_SMG(PC pc, PetscOptionItems PetscOptionsObject)
3220: {
3221:   PC_SMG *ex = (PC_SMG *)pc->data;

3223:   PetscFunctionBegin;
3224:   PetscOptionsHeadBegin(PetscOptionsObject, "SMG options");

3226:   PetscCall(PetscOptionsInt("-pc_smg_its", "Number of iterations of SMG to use as preconditioner", "HYPRE_StructSMGSetMaxIter", ex->its, &ex->its, NULL));
3227:   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));
3228:   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));
3229:   PetscCall(PetscOptionsReal("-pc_smg_tol", "Tolerance of SMG", "HYPRE_StructSMGSetTol", ex->tol, &ex->tol, NULL));

3231:   PetscOptionsHeadEnd();
3232:   PetscFunctionReturn(PETSC_SUCCESS);
3233: }

3235: static PetscErrorCode PCApply_SMG(PC pc, Vec x, Vec y)
3236: {
3237:   PC_SMG            *ex = (PC_SMG *)pc->data;
3238:   PetscScalar       *yy;
3239:   const PetscScalar *xx;
3240:   PetscInt           ilower[3], iupper[3];
3241:   HYPRE_Int          hlower[3], hupper[3];
3242:   Mat_HYPREStruct   *mx = (Mat_HYPREStruct *)pc->pmat->data;

3244:   PetscFunctionBegin;
3245:   PetscCall(PetscCitationsRegister(hypreCitation, &cite));
3246:   PetscCall(DMDAGetCorners(mx->da, &ilower[0], &ilower[1], &ilower[2], &iupper[0], &iupper[1], &iupper[2]));
3247:   /* when HYPRE_MIXEDINT is defined, sizeof(HYPRE_Int) == 32 */
3248:   iupper[0] += ilower[0] - 1;
3249:   iupper[1] += ilower[1] - 1;
3250:   iupper[2] += ilower[2] - 1;
3251:   hlower[0] = (HYPRE_Int)ilower[0];
3252:   hlower[1] = (HYPRE_Int)ilower[1];
3253:   hlower[2] = (HYPRE_Int)ilower[2];
3254:   hupper[0] = (HYPRE_Int)iupper[0];
3255:   hupper[1] = (HYPRE_Int)iupper[1];
3256:   hupper[2] = (HYPRE_Int)iupper[2];

3258:   /* copy x values over to hypre */
3259:   PetscCallExternal(HYPRE_StructVectorSetConstantValues, mx->hb, 0.0);
3260:   PetscCall(VecGetArrayRead(x, &xx));
3261:   PetscCallExternal(HYPRE_StructVectorSetBoxValues, mx->hb, hlower, hupper, (HYPRE_Complex *)xx);
3262:   PetscCall(VecRestoreArrayRead(x, &xx));
3263:   PetscCallExternal(HYPRE_StructVectorAssemble, mx->hb);
3264:   PetscCallExternal(HYPRE_StructSMGSolve, ex->hsolver, mx->hmat, mx->hb, mx->hx);

3266:   /* copy solution values back to PETSc */
3267:   PetscCall(VecGetArray(y, &yy));
3268:   PetscCallExternal(HYPRE_StructVectorGetBoxValues, mx->hx, hlower, hupper, (HYPRE_Complex *)yy);
3269:   PetscCall(VecRestoreArray(y, &yy));
3270:   PetscFunctionReturn(PETSC_SUCCESS);
3271: }

3273: 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)
3274: {
3275:   PC_SMG   *jac = (PC_SMG *)pc->data;
3276:   HYPRE_Int oits;

3278:   PetscFunctionBegin;
3279:   PetscCall(PetscCitationsRegister(hypreCitation, &cite));
3280:   PetscCallExternal(HYPRE_StructSMGSetMaxIter, jac->hsolver, its * jac->its);
3281:   PetscCallExternal(HYPRE_StructSMGSetTol, jac->hsolver, rtol);

3283:   PetscCall(PCApply_SMG(pc, b, y));
3284:   PetscCallExternal(HYPRE_StructSMGGetNumIterations, jac->hsolver, &oits);
3285:   *outits = oits;
3286:   if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS;
3287:   else *reason = PCRICHARDSON_CONVERGED_RTOL;
3288:   PetscCallExternal(HYPRE_StructSMGSetTol, jac->hsolver, jac->tol);
3289:   PetscCallExternal(HYPRE_StructSMGSetMaxIter, jac->hsolver, jac->its);
3290:   PetscFunctionReturn(PETSC_SUCCESS);
3291: }

3293: static PetscErrorCode PCSetUp_SMG(PC pc)
3294: {
3295:   PetscInt         i, dim;
3296:   PC_SMG          *ex = (PC_SMG *)pc->data;
3297:   Mat_HYPREStruct *mx = (Mat_HYPREStruct *)pc->pmat->data;
3298:   PetscBool        flg;
3299:   DMBoundaryType   p[3];
3300:   PetscInt         M[3];

3302:   PetscFunctionBegin;
3303:   PetscCall(PetscObjectTypeCompare((PetscObject)pc->pmat, MATHYPRESTRUCT, &flg));
3304:   PetscCheck(flg, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_INCOMP, "Must use MATHYPRESTRUCT with this preconditioner");

3306:   PetscCall(DMDAGetInfo(mx->da, &dim, &M[0], &M[1], &M[2], 0, 0, 0, 0, 0, &p[0], &p[1], &p[2], 0));
3307:   // Check if power of 2 in periodic directions
3308:   for (i = 0; i < dim; i++) {
3309:     PetscCheck((M[i] & (M[i] - 1)) == 0 || p[i] != DM_BOUNDARY_PERIODIC, 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]);
3310:   }

3312:   /* create the hypre solver object and set its information */
3313:   if (ex->hsolver) PetscCallExternal(HYPRE_StructSMGDestroy, ex->hsolver);
3314:   PetscCallExternal(HYPRE_StructSMGCreate, ex->hcomm, &ex->hsolver);
3315:   // The hypre options must be set here and not in SetFromOptions because it is created here!
3316:   PetscCallExternal(HYPRE_StructSMGSetMaxIter, ex->hsolver, ex->its);
3317:   PetscCallExternal(HYPRE_StructSMGSetNumPreRelax, ex->hsolver, ex->num_pre_relax);
3318:   PetscCallExternal(HYPRE_StructSMGSetNumPostRelax, ex->hsolver, ex->num_post_relax);
3319:   PetscCallExternal(HYPRE_StructSMGSetTol, ex->hsolver, ex->tol);

3321:   PetscCallExternal(HYPRE_StructSMGSetup, ex->hsolver, mx->hmat, mx->hb, mx->hx);
3322:   PetscCallExternal(HYPRE_StructSMGSetZeroGuess, ex->hsolver);
3323:   PetscFunctionReturn(PETSC_SUCCESS);
3324: }

3326: /*MC
3327:      PCSMG - the hypre (structured grid) SMG multigrid solver

3329:    Level: advanced

3331:    Options Database Keys:
3332: + -pc_smg_its <its> - number of iterations of SMG to use as preconditioner
3333: . -pc_smg_num_pre_relax <steps> - number of smoothing steps before coarse grid
3334: . -pc_smg_num_post_relax <steps> - number of smoothing steps after coarse grid
3335: - -pc_smg_tol <tol> - tolerance of SMG

3337:    Notes:
3338:    This is for CELL-centered descretizations

3340:    This must be used with the `MATHYPRESTRUCT` `MatType`.

3342:    This does not provide all the functionality of  hypre's SMG solver, it supports only one block per process defined by a PETSc `DMDA`.

3344:    See `PCSYSPFMG`, `PCSMG`, `PCPFMG`, and `PCHYPRE` for access to hypre's other preconditioners

3346: .seealso:  `PCMG`, `MATHYPRESTRUCT`, `PCPFMG`, `PCSYSPFMG`, `PCHYPRE`, `PCGAMG`
3347: M*/

3349: PETSC_EXTERN PetscErrorCode PCCreate_SMG(PC pc)
3350: {
3351:   PC_SMG *ex;

3353:   PetscFunctionBegin;
3354:   PetscCall(PetscNew(&ex));
3355:   pc->data = ex;

3357:   ex->its            = 1;
3358:   ex->tol            = 1.e-8;
3359:   ex->num_pre_relax  = 1;
3360:   ex->num_post_relax = 1;

3362:   pc->ops->setfromoptions  = PCSetFromOptions_SMG;
3363:   pc->ops->view            = PCView_SMG;
3364:   pc->ops->destroy         = PCDestroy_SMG;
3365:   pc->ops->apply           = PCApply_SMG;
3366:   pc->ops->applyrichardson = PCApplyRichardson_SMG;
3367:   pc->ops->setup           = PCSetUp_SMG;

3369:   PetscCall(PetscCommGetComm(PetscObjectComm((PetscObject)pc), &ex->hcomm));
3370:   PetscHYPREInitialize();
3371:   PetscCallExternal(HYPRE_StructSMGCreate, ex->hcomm, &ex->hsolver);
3372:   PetscFunctionReturn(PETSC_SUCCESS);
3373: }