xref: /petsc/src/ksp/pc/impls/hmg/hmg.c (revision 48a46eb9bd028bec07ec0f396b1a3abb43f14558)
1360ee056SFande Kong #include <petscdm.h>
2360ee056SFande Kong #include <petscctable.h>
3360ee056SFande Kong #include <petsc/private/matimpl.h>
4360ee056SFande Kong #include <petsc/private/pcmgimpl.h>
5360ee056SFande Kong #include <petsc/private/pcimpl.h> /*I "petscpc.h" I*/
6360ee056SFande Kong 
7360ee056SFande Kong typedef struct {
88a2c336bSFande Kong   PC        innerpc;       /* A MG inner PC (Hypre or PCGAMG) to setup interpolations and coarse operators  */
98a2c336bSFande Kong   char     *innerpctype;   /* PCGAMG or PCHYPRE */
108a2c336bSFande Kong   PetscBool reuseinterp;   /* A flag indicates if or not to reuse the interpolations */
1149c604d5SFande Kong   PetscBool subcoarsening; /* If or not to use a subspace-based coarsening algorithm */
1249c604d5SFande Kong   PetscBool usematmaij;    /* If or not to use MatMAIJ for saving memory */
1349c604d5SFande Kong   PetscInt  component;     /* Which subspace is used for the subspace-based coarsening algorithm? */
14360ee056SFande Kong } PC_HMG;
15360ee056SFande Kong 
16dbbe0bcdSBarry Smith PetscErrorCode PCSetFromOptions_HMG(PC, PetscOptionItems *);
17360ee056SFande Kong PetscErrorCode PCReset_MG(PC);
188a2c336bSFande Kong 
199371c9d4SSatish Balay static PetscErrorCode PCHMGExtractSubMatrix_Private(Mat pmat, Mat *submat, MatReuse reuse, PetscInt component, PetscInt blocksize) {
20360ee056SFande Kong   IS       isrow;
218a2c336bSFande Kong   PetscInt rstart, rend;
22360ee056SFande Kong   MPI_Comm comm;
23360ee056SFande Kong 
24360ee056SFande Kong   PetscFunctionBegin;
259566063dSJacob Faibussowitsch   PetscCall(PetscObjectGetComm((PetscObject)pmat, &comm));
2663a3b9bcSJacob Faibussowitsch   PetscCheck(component < blocksize, comm, PETSC_ERR_ARG_INCOMP, "Component %" PetscInt_FMT " should be less than block size %" PetscInt_FMT " ", component, blocksize);
279566063dSJacob Faibussowitsch   PetscCall(MatGetOwnershipRange(pmat, &rstart, &rend));
2863a3b9bcSJacob Faibussowitsch   PetscCheck((rend - rstart) % blocksize == 0, comm, PETSC_ERR_ARG_INCOMP, "Block size %" PetscInt_FMT " is inconsistent for [%" PetscInt_FMT ", %" PetscInt_FMT ") ", blocksize, rstart, rend);
299566063dSJacob Faibussowitsch   PetscCall(ISCreateStride(comm, (rend - rstart) / blocksize, rstart + component, blocksize, &isrow));
309566063dSJacob Faibussowitsch   PetscCall(MatCreateSubMatrix(pmat, isrow, isrow, reuse, submat));
319566063dSJacob Faibussowitsch   PetscCall(ISDestroy(&isrow));
32360ee056SFande Kong   PetscFunctionReturn(0);
33360ee056SFande Kong }
34360ee056SFande Kong 
359371c9d4SSatish Balay static PetscErrorCode PCHMGExpandInterpolation_Private(Mat subinterp, Mat *interp, PetscInt blocksize) {
36360ee056SFande Kong   PetscInt           subrstart, subrend, subrowsize, subcolsize, subcstart, subcend, rowsize, colsize;
378a2c336bSFande Kong   PetscInt           subrow, row, nz, *d_nnz, *o_nnz, i, j, dnz, onz, max_nz, *indices;
388a2c336bSFande Kong   const PetscInt    *idx;
398a2c336bSFande Kong   const PetscScalar *values;
408a2c336bSFande Kong   MPI_Comm           comm;
41360ee056SFande Kong 
42360ee056SFande Kong   PetscFunctionBegin;
439566063dSJacob Faibussowitsch   PetscCall(PetscObjectGetComm((PetscObject)subinterp, &comm));
449566063dSJacob Faibussowitsch   PetscCall(MatGetOwnershipRange(subinterp, &subrstart, &subrend));
458a2c336bSFande Kong   subrowsize = subrend - subrstart;
46360ee056SFande Kong   rowsize    = subrowsize * blocksize;
479566063dSJacob Faibussowitsch   PetscCall(PetscCalloc2(rowsize, &d_nnz, rowsize, &o_nnz));
489566063dSJacob Faibussowitsch   PetscCall(MatGetOwnershipRangeColumn(subinterp, &subcstart, &subcend));
498a2c336bSFande Kong   subcolsize = subcend - subcstart;
508a2c336bSFande Kong   colsize    = subcolsize * blocksize;
51360ee056SFande Kong   max_nz     = 0;
528a2c336bSFande Kong   for (subrow = subrstart; subrow < subrend; subrow++) {
539566063dSJacob Faibussowitsch     PetscCall(MatGetRow(subinterp, subrow, &nz, &idx, NULL));
54360ee056SFande Kong     if (max_nz < nz) max_nz = nz;
559371c9d4SSatish Balay     dnz = 0;
569371c9d4SSatish Balay     onz = 0;
57360ee056SFande Kong     for (i = 0; i < nz; i++) {
588a2c336bSFande Kong       if (idx[i] >= subcstart && idx[i] < subcend) dnz++;
598a2c336bSFande Kong       else onz++;
60360ee056SFande Kong     }
61360ee056SFande Kong     for (i = 0; i < blocksize; i++) {
62360ee056SFande Kong       d_nnz[(subrow - subrstart) * blocksize + i] = dnz;
63360ee056SFande Kong       o_nnz[(subrow - subrstart) * blocksize + i] = onz;
64360ee056SFande Kong     }
659566063dSJacob Faibussowitsch     PetscCall(MatRestoreRow(subinterp, subrow, &nz, &idx, NULL));
66360ee056SFande Kong   }
679566063dSJacob Faibussowitsch   PetscCall(MatCreateAIJ(comm, rowsize, colsize, PETSC_DETERMINE, PETSC_DETERMINE, 0, d_nnz, 0, o_nnz, interp));
689566063dSJacob Faibussowitsch   PetscCall(MatSetOption(*interp, MAT_IGNORE_OFF_PROC_ENTRIES, PETSC_TRUE));
699566063dSJacob Faibussowitsch   PetscCall(MatSetOption(*interp, MAT_IGNORE_ZERO_ENTRIES, PETSC_TRUE));
709566063dSJacob Faibussowitsch   PetscCall(MatSetOption(*interp, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_TRUE));
719566063dSJacob Faibussowitsch   PetscCall(MatSetFromOptions(*interp));
72360ee056SFande Kong 
739566063dSJacob Faibussowitsch   PetscCall(MatSetUp(*interp));
749566063dSJacob Faibussowitsch   PetscCall(PetscFree2(d_nnz, o_nnz));
759566063dSJacob Faibussowitsch   PetscCall(PetscMalloc1(max_nz, &indices));
768a2c336bSFande Kong   for (subrow = subrstart; subrow < subrend; subrow++) {
779566063dSJacob Faibussowitsch     PetscCall(MatGetRow(subinterp, subrow, &nz, &idx, &values));
78360ee056SFande Kong     for (i = 0; i < blocksize; i++) {
79360ee056SFande Kong       row = subrow * blocksize + i;
809371c9d4SSatish Balay       for (j = 0; j < nz; j++) { indices[j] = idx[j] * blocksize + i; }
819566063dSJacob Faibussowitsch       PetscCall(MatSetValues(*interp, 1, &row, nz, indices, values, INSERT_VALUES));
82360ee056SFande Kong     }
839566063dSJacob Faibussowitsch     PetscCall(MatRestoreRow(subinterp, subrow, &nz, &idx, &values));
84360ee056SFande Kong   }
859566063dSJacob Faibussowitsch   PetscCall(PetscFree(indices));
869566063dSJacob Faibussowitsch   PetscCall(MatAssemblyBegin(*interp, MAT_FINAL_ASSEMBLY));
879566063dSJacob Faibussowitsch   PetscCall(MatAssemblyEnd(*interp, MAT_FINAL_ASSEMBLY));
88360ee056SFande Kong   PetscFunctionReturn(0);
89360ee056SFande Kong }
90360ee056SFande Kong 
919371c9d4SSatish Balay PetscErrorCode PCSetUp_HMG(PC pc) {
92360ee056SFande Kong   Mat              PA, submat;
93360ee056SFande Kong   PC_MG           *mg  = (PC_MG *)pc->data;
94360ee056SFande Kong   PC_HMG          *hmg = (PC_HMG *)mg->innerctx;
95360ee056SFande Kong   MPI_Comm         comm;
96360ee056SFande Kong   PetscInt         level;
97360ee056SFande Kong   PetscInt         num_levels;
988a2c336bSFande Kong   Mat             *operators, *interpolations;
998a2c336bSFande Kong   PetscInt         blocksize;
100fd2dd295SFande Kong   const char      *prefix;
10107a4832bSFande Kong   PCMGGalerkinType galerkin;
102360ee056SFande Kong 
103360ee056SFande Kong   PetscFunctionBegin;
1049566063dSJacob Faibussowitsch   PetscCall(PetscObjectGetComm((PetscObject)pc, &comm));
105360ee056SFande Kong   if (pc->setupcalled) {
1065f36d8f2SFande Kong     if (hmg->reuseinterp) {
1075f36d8f2SFande Kong       /* If we did not use Galerkin in the last call or we have a different sparsity pattern now,
1085f36d8f2SFande Kong       * we have to build from scratch
10907a4832bSFande Kong       * */
1109566063dSJacob Faibussowitsch       PetscCall(PCMGGetGalerkin(pc, &galerkin));
1115f36d8f2SFande Kong       if (galerkin == PC_MG_GALERKIN_NONE || pc->flag != SAME_NONZERO_PATTERN) pc->setupcalled = PETSC_FALSE;
1129566063dSJacob Faibussowitsch       PetscCall(PCMGSetGalerkin(pc, PC_MG_GALERKIN_PMAT));
1139566063dSJacob Faibussowitsch       PetscCall(PCSetUp_MG(pc));
114360ee056SFande Kong       PetscFunctionReturn(0);
115360ee056SFande Kong     } else {
1169566063dSJacob Faibussowitsch       PetscCall(PCReset_MG(pc));
117360ee056SFande Kong       pc->setupcalled = PETSC_FALSE;
118360ee056SFande Kong     }
119360ee056SFande Kong   }
120360ee056SFande Kong 
1218a2c336bSFande Kong   /* Create an inner PC (GAMG or HYPRE) */
1228a2c336bSFande Kong   if (!hmg->innerpc) {
1239566063dSJacob Faibussowitsch     PetscCall(PCCreate(comm, &hmg->innerpc));
124fd2dd295SFande Kong     /* If users do not set an inner pc type, we need to set a default value */
125fd2dd295SFande Kong     if (!hmg->innerpctype) {
126fd2dd295SFande Kong       /* If hypre is available, use hypre, otherwise, use gamg */
127fd2dd295SFande Kong #if PETSC_HAVE_HYPRE
1289566063dSJacob Faibussowitsch       PetscCall(PetscStrallocpy(PCHYPRE, &(hmg->innerpctype)));
129fd2dd295SFande Kong #else
1309566063dSJacob Faibussowitsch       PetscCall(PetscStrallocpy(PCGAMG, &(hmg->innerpctype)));
131fd2dd295SFande Kong #endif
1328a2c336bSFande Kong     }
1339566063dSJacob Faibussowitsch     PetscCall(PCSetType(hmg->innerpc, hmg->innerpctype));
134360ee056SFande Kong   }
1359566063dSJacob Faibussowitsch   PetscCall(PCGetOperators(pc, NULL, &PA));
1368a2c336bSFande Kong   /* Users need to correctly set a block size of matrix in order to use subspace coarsening */
1379566063dSJacob Faibussowitsch   PetscCall(MatGetBlockSize(PA, &blocksize));
1388a2c336bSFande Kong   if (blocksize <= 1) hmg->subcoarsening = PETSC_FALSE;
1398a2c336bSFande Kong   /* Extract a submatrix for constructing subinterpolations */
1408a2c336bSFande Kong   if (hmg->subcoarsening) {
1419566063dSJacob Faibussowitsch     PetscCall(PCHMGExtractSubMatrix_Private(PA, &submat, MAT_INITIAL_MATRIX, hmg->component, blocksize));
142360ee056SFande Kong     PA = submat;
143360ee056SFande Kong   }
1449566063dSJacob Faibussowitsch   PetscCall(PCSetOperators(hmg->innerpc, PA, PA));
145*48a46eb9SPierre Jolivet   if (hmg->subcoarsening) PetscCall(MatDestroy(&PA));
1468a2c336bSFande Kong   /* Setup inner PC correctly. During this step, matrix will be coarsened */
1479566063dSJacob Faibussowitsch   PetscCall(PCSetUseAmat(hmg->innerpc, PETSC_FALSE));
1489566063dSJacob Faibussowitsch   PetscCall(PetscObjectGetOptionsPrefix((PetscObject)pc, &prefix));
1499566063dSJacob Faibussowitsch   PetscCall(PetscObjectSetOptionsPrefix((PetscObject)hmg->innerpc, prefix));
1509566063dSJacob Faibussowitsch   PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)hmg->innerpc, "hmg_inner_"));
1519566063dSJacob Faibussowitsch   PetscCall(PCSetFromOptions(hmg->innerpc));
1529566063dSJacob Faibussowitsch   PetscCall(PCSetUp(hmg->innerpc));
153360ee056SFande Kong 
1548a2c336bSFande Kong   /* Obtain interpolations IN PLACE. For BoomerAMG, (I,J,data) is reused to avoid memory overhead */
1559566063dSJacob Faibussowitsch   PetscCall(PCGetInterpolations(hmg->innerpc, &num_levels, &interpolations));
1568a2c336bSFande Kong   /* We can reuse the coarse operators when we do the full space coarsening */
157*48a46eb9SPierre Jolivet   if (!hmg->subcoarsening) PetscCall(PCGetCoarseOperators(hmg->innerpc, &num_levels, &operators));
158360ee056SFande Kong 
1599566063dSJacob Faibussowitsch   PetscCall(PCDestroy(&hmg->innerpc));
1600a545947SLisandro Dalcin   hmg->innerpc = NULL;
1619566063dSJacob Faibussowitsch   PetscCall(PCMGSetLevels_MG(pc, num_levels, NULL));
1628a2c336bSFande Kong   /* Set coarse matrices and interpolations to PCMG */
163360ee056SFande Kong   for (level = num_levels - 1; level > 0; level--) {
1640a545947SLisandro Dalcin     Mat P = NULL, pmat = NULL;
165360ee056SFande Kong     Vec b, x, r;
1668a2c336bSFande Kong     if (hmg->subcoarsening) {
1674bb91820SFande Kong       if (hmg->usematmaij) {
1689566063dSJacob Faibussowitsch         PetscCall(MatCreateMAIJ(interpolations[level - 1], blocksize, &P));
1699566063dSJacob Faibussowitsch         PetscCall(MatDestroy(&interpolations[level - 1]));
1704bb91820SFande Kong       } else {
1718a2c336bSFande Kong         /* Grow interpolation. In the future, we should use MAIJ */
1729566063dSJacob Faibussowitsch         PetscCall(PCHMGExpandInterpolation_Private(interpolations[level - 1], &P, blocksize));
1739566063dSJacob Faibussowitsch         PetscCall(MatDestroy(&interpolations[level - 1]));
1744bb91820SFande Kong       }
175360ee056SFande Kong     } else {
1768a2c336bSFande Kong       P = interpolations[level - 1];
177360ee056SFande Kong     }
1789566063dSJacob Faibussowitsch     PetscCall(MatCreateVecs(P, &b, &r));
1799566063dSJacob Faibussowitsch     PetscCall(PCMGSetInterpolation(pc, level, P));
1809566063dSJacob Faibussowitsch     PetscCall(PCMGSetRestriction(pc, level, P));
1819566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&P));
1828a2c336bSFande Kong     /* We reuse the matrices when we do not do subspace coarsening */
1838a2c336bSFande Kong     if ((level - 1) >= 0 && !hmg->subcoarsening) {
1848a2c336bSFande Kong       pmat = operators[level - 1];
1859566063dSJacob Faibussowitsch       PetscCall(PCMGSetOperators(pc, level - 1, pmat, pmat));
1869566063dSJacob Faibussowitsch       PetscCall(MatDestroy(&pmat));
187360ee056SFande Kong     }
1889566063dSJacob Faibussowitsch     PetscCall(PCMGSetRhs(pc, level - 1, b));
189360ee056SFande Kong 
1909566063dSJacob Faibussowitsch     PetscCall(PCMGSetR(pc, level, r));
1919566063dSJacob Faibussowitsch     PetscCall(VecDestroy(&r));
192360ee056SFande Kong 
1939566063dSJacob Faibussowitsch     PetscCall(VecDuplicate(b, &x));
1949566063dSJacob Faibussowitsch     PetscCall(PCMGSetX(pc, level - 1, x));
1959566063dSJacob Faibussowitsch     PetscCall(VecDestroy(&x));
1969566063dSJacob Faibussowitsch     PetscCall(VecDestroy(&b));
197360ee056SFande Kong   }
1989566063dSJacob Faibussowitsch   PetscCall(PetscFree(interpolations));
199*48a46eb9SPierre Jolivet   if (!hmg->subcoarsening) PetscCall(PetscFree(operators));
2008a2c336bSFande Kong   /* Turn Galerkin off when we already have coarse operators */
2019566063dSJacob Faibussowitsch   PetscCall(PCMGSetGalerkin(pc, hmg->subcoarsening ? PC_MG_GALERKIN_PMAT : PC_MG_GALERKIN_NONE));
2029566063dSJacob Faibussowitsch   PetscCall(PCSetDM(pc, NULL));
2039566063dSJacob Faibussowitsch   PetscCall(PCSetUseAmat(pc, PETSC_FALSE));
204d0609cedSBarry Smith   PetscObjectOptionsBegin((PetscObject)pc);
205dbbe0bcdSBarry Smith   PetscCall(PCSetFromOptions_MG(pc, PetscOptionsObject)); /* should be called in PCSetFromOptions_HMG(), but cannot be called prior to PCMGSetLevels() */
206d0609cedSBarry Smith   PetscOptionsEnd();
2079566063dSJacob Faibussowitsch   PetscCall(PCSetUp_MG(pc));
208360ee056SFande Kong   PetscFunctionReturn(0);
209360ee056SFande Kong }
210360ee056SFande Kong 
2119371c9d4SSatish Balay PetscErrorCode PCDestroy_HMG(PC pc) {
212360ee056SFande Kong   PC_MG  *mg  = (PC_MG *)pc->data;
2138a2c336bSFande Kong   PC_HMG *hmg = (PC_HMG *)mg->innerctx;
214360ee056SFande Kong 
215360ee056SFande Kong   PetscFunctionBegin;
2169566063dSJacob Faibussowitsch   PetscCall(PCDestroy(&hmg->innerpc));
2179566063dSJacob Faibussowitsch   PetscCall(PetscFree(hmg->innerpctype));
2189566063dSJacob Faibussowitsch   PetscCall(PetscFree(hmg));
2199566063dSJacob Faibussowitsch   PetscCall(PCDestroy_MG(pc));
220fd2dd295SFande Kong 
2219566063dSJacob Faibussowitsch   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHMGSetReuseInterpolation_C", NULL));
2229566063dSJacob Faibussowitsch   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHMGSetUseSubspaceCoarsening_C", NULL));
2239566063dSJacob Faibussowitsch   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHMGSetInnerPCType_C", NULL));
2249566063dSJacob Faibussowitsch   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHMGSetCoarseningComponent_C", NULL));
2252e956fe4SStefano Zampini   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHMGUseMatMAIJ_C", NULL));
226360ee056SFande Kong   PetscFunctionReturn(0);
227360ee056SFande Kong }
228360ee056SFande Kong 
2299371c9d4SSatish Balay PetscErrorCode PCView_HMG(PC pc, PetscViewer viewer) {
230360ee056SFande Kong   PC_MG    *mg  = (PC_MG *)pc->data;
2318a2c336bSFande Kong   PC_HMG   *hmg = (PC_HMG *)mg->innerctx;
232360ee056SFande Kong   PetscBool iascii;
233360ee056SFande Kong 
234360ee056SFande Kong   PetscFunctionBegin;
2359566063dSJacob Faibussowitsch   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
236360ee056SFande Kong   if (iascii) {
2379566063dSJacob Faibussowitsch     PetscCall(PetscViewerASCIIPrintf(viewer, " Reuse interpolation: %s\n", hmg->reuseinterp ? "true" : "false"));
2389566063dSJacob Faibussowitsch     PetscCall(PetscViewerASCIIPrintf(viewer, " Use subspace coarsening: %s\n", hmg->subcoarsening ? "true" : "false"));
23963a3b9bcSJacob Faibussowitsch     PetscCall(PetscViewerASCIIPrintf(viewer, " Coarsening component: %" PetscInt_FMT " \n", hmg->component));
2409566063dSJacob Faibussowitsch     PetscCall(PetscViewerASCIIPrintf(viewer, " Use MatMAIJ: %s \n", hmg->usematmaij ? "true" : "false"));
2419566063dSJacob Faibussowitsch     PetscCall(PetscViewerASCIIPrintf(viewer, " Inner PC type: %s \n", hmg->innerpctype));
242360ee056SFande Kong   }
2439566063dSJacob Faibussowitsch   PetscCall(PCView_MG(pc, viewer));
244360ee056SFande Kong   PetscFunctionReturn(0);
245360ee056SFande Kong }
246360ee056SFande Kong 
2479371c9d4SSatish Balay PetscErrorCode PCSetFromOptions_HMG(PC pc, PetscOptionItems *PetscOptionsObject) {
248360ee056SFande Kong   PC_MG  *mg  = (PC_MG *)pc->data;
2498a2c336bSFande Kong   PC_HMG *hmg = (PC_HMG *)mg->innerctx;
250360ee056SFande Kong 
251360ee056SFande Kong   PetscFunctionBegin;
252d0609cedSBarry Smith   PetscOptionsHeadBegin(PetscOptionsObject, "HMG");
2539566063dSJacob Faibussowitsch   PetscCall(PetscOptionsBool("-pc_hmg_reuse_interpolation", "Reuse the interpolation operators when possible (cheaper, weaker when matrix entries change a lot)", "PCHMGSetReuseInterpolation", hmg->reuseinterp, &hmg->reuseinterp, NULL));
2549566063dSJacob Faibussowitsch   PetscCall(PetscOptionsBool("-pc_hmg_use_subspace_coarsening", "Use the subspace coarsening to compute the interpolations", "PCHMGSetUseSubspaceCoarsening", hmg->subcoarsening, &hmg->subcoarsening, NULL));
2559566063dSJacob Faibussowitsch   PetscCall(PetscOptionsBool("-pc_hmg_use_matmaij", "Use MatMAIJ store interpolation for saving memory", "PCHMGSetInnerPCType", hmg->usematmaij, &hmg->usematmaij, NULL));
2569566063dSJacob Faibussowitsch   PetscCall(PetscOptionsInt("-pc_hmg_coarsening_component", "Which component is chosen for the subspace-based coarsening algorithm", "PCHMGSetCoarseningComponent", hmg->component, &hmg->component, NULL));
257d0609cedSBarry Smith   PetscOptionsHeadEnd();
258360ee056SFande Kong   PetscFunctionReturn(0);
259360ee056SFande Kong }
260360ee056SFande Kong 
2619371c9d4SSatish Balay static PetscErrorCode PCHMGSetReuseInterpolation_HMG(PC pc, PetscBool reuse) {
26207a4832bSFande Kong   PC_MG  *mg  = (PC_MG *)pc->data;
26307a4832bSFande Kong   PC_HMG *hmg = (PC_HMG *)mg->innerctx;
264fd2dd295SFande Kong 
265fd2dd295SFande Kong   PetscFunctionBegin;
266fd2dd295SFande Kong   hmg->reuseinterp = reuse;
267fd2dd295SFande Kong   PetscFunctionReturn(0);
268fd2dd295SFande Kong }
269fd2dd295SFande Kong 
270b155ba7fSFande Kong /*@
2718a2c336bSFande Kong    PCHMGSetReuseInterpolation - Reuse interpolation matrices in HMG
2728a2c336bSFande Kong 
2738a2c336bSFande Kong    Logically Collective on PC
2748a2c336bSFande Kong 
2758a2c336bSFande Kong    Input Parameters:
2768a2c336bSFande Kong +  pc - the HMG context
2778a2c336bSFande Kong -  reuse - True indicates that HMG will reuse the interpolations
2788a2c336bSFande Kong 
2798a2c336bSFande Kong    Options Database Keys:
28049c604d5SFande Kong .  -pc_hmg_reuse_interpolation <true | false> - Whether or not to reuse the interpolations. If true, it potentially save the compute time.
2818a2c336bSFande Kong 
2828a2c336bSFande Kong    Level: beginner
2838a2c336bSFande Kong 
2848a2c336bSFande Kong .keywords: HMG, multigrid, interpolation, reuse, set
2858a2c336bSFande Kong 
286db781477SPatrick Sanan .seealso: `PCHMG`
287b155ba7fSFande Kong @*/
2889371c9d4SSatish Balay PetscErrorCode PCHMGSetReuseInterpolation(PC pc, PetscBool reuse) {
289fd2dd295SFande Kong   PetscFunctionBegin;
290fd2dd295SFande Kong   PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
291cac4c232SBarry Smith   PetscUseMethod(pc, "PCHMGSetReuseInterpolation_C", (PC, PetscBool), (pc, reuse));
292fd2dd295SFande Kong   PetscFunctionReturn(0);
293fd2dd295SFande Kong }
294fd2dd295SFande Kong 
2959371c9d4SSatish Balay static PetscErrorCode PCHMGSetUseSubspaceCoarsening_HMG(PC pc, PetscBool subspace) {
29607a4832bSFande Kong   PC_MG  *mg  = (PC_MG *)pc->data;
29707a4832bSFande Kong   PC_HMG *hmg = (PC_HMG *)mg->innerctx;
2988a2c336bSFande Kong 
2998a2c336bSFande Kong   PetscFunctionBegin;
300fd2dd295SFande Kong   hmg->subcoarsening = subspace;
3018a2c336bSFande Kong   PetscFunctionReturn(0);
3028a2c336bSFande Kong }
3038a2c336bSFande Kong 
304b155ba7fSFande Kong /*@
3058a2c336bSFande Kong    PCHMGSetUseSubspaceCoarsening - Use subspace coarsening in HMG
3068a2c336bSFande Kong 
3078a2c336bSFande Kong    Logically Collective on PC
3088a2c336bSFande Kong 
3098a2c336bSFande Kong    Input Parameters:
3108a2c336bSFande Kong +  pc - the HMG context
3118a2c336bSFande Kong -  reuse - True indicates that HMG will use the subspace coarsening
3128a2c336bSFande Kong 
3138a2c336bSFande Kong    Options Database Keys:
31449c604d5SFande Kong .  -pc_hmg_use_subspace_coarsening  <true | false> - Whether or not to use subspace coarsening (that is, coarsen a submatrix).
3158a2c336bSFande Kong 
3168a2c336bSFande Kong    Level: beginner
3178a2c336bSFande Kong 
3188a2c336bSFande Kong .keywords: HMG, multigrid, interpolation, subspace, coarsening
3198a2c336bSFande Kong 
320db781477SPatrick Sanan .seealso: `PCHMG`
321b155ba7fSFande Kong @*/
3229371c9d4SSatish Balay PetscErrorCode PCHMGSetUseSubspaceCoarsening(PC pc, PetscBool subspace) {
3238a2c336bSFande Kong   PetscFunctionBegin;
3248a2c336bSFande Kong   PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
325cac4c232SBarry Smith   PetscUseMethod(pc, "PCHMGSetUseSubspaceCoarsening_C", (PC, PetscBool), (pc, subspace));
326fd2dd295SFande Kong   PetscFunctionReturn(0);
327fd2dd295SFande Kong }
328fd2dd295SFande Kong 
3299371c9d4SSatish Balay static PetscErrorCode PCHMGSetInnerPCType_HMG(PC pc, PCType type) {
33007a4832bSFande Kong   PC_MG  *mg  = (PC_MG *)pc->data;
33107a4832bSFande Kong   PC_HMG *hmg = (PC_HMG *)mg->innerctx;
332fd2dd295SFande Kong 
333fd2dd295SFande Kong   PetscFunctionBegin;
3349566063dSJacob Faibussowitsch   PetscCall(PetscStrallocpy(type, &(hmg->innerpctype)));
3358a2c336bSFande Kong   PetscFunctionReturn(0);
3368a2c336bSFande Kong }
3378a2c336bSFande Kong 
338b155ba7fSFande Kong /*@C
3398a2c336bSFande Kong    PCHMGSetInnerPCType - Set an inner PC type
3408a2c336bSFande Kong 
3418a2c336bSFande Kong    Logically Collective on PC
3428a2c336bSFande Kong 
3438a2c336bSFande Kong    Input Parameters:
3448a2c336bSFande Kong +  pc - the HMG context
3458a2c336bSFande Kong -  type - <hypre, gamg> coarsening algorithm
3468a2c336bSFande Kong 
3478a2c336bSFande Kong    Options Database Keys:
34849c604d5SFande Kong .  -hmg_inner_pc_type <hypre, gamg> - What method is used to coarsen matrix
3498a2c336bSFande Kong 
3508a2c336bSFande Kong    Level: beginner
3518a2c336bSFande Kong 
3528a2c336bSFande Kong .keywords: HMG, multigrid, interpolation, coarsening
3538a2c336bSFande Kong 
354db781477SPatrick Sanan .seealso: `PCHMG`, `PCType`
355b155ba7fSFande Kong @*/
3569371c9d4SSatish Balay PetscErrorCode PCHMGSetInnerPCType(PC pc, PCType type) {
3578a2c336bSFande Kong   PetscFunctionBegin;
3588a2c336bSFande Kong   PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
359cac4c232SBarry Smith   PetscUseMethod(pc, "PCHMGSetInnerPCType_C", (PC, PCType), (pc, type));
3608a2c336bSFande Kong   PetscFunctionReturn(0);
3618a2c336bSFande Kong }
3628a2c336bSFande Kong 
3639371c9d4SSatish Balay static PetscErrorCode PCHMGSetCoarseningComponent_HMG(PC pc, PetscInt component) {
36449c604d5SFande Kong   PC_MG  *mg  = (PC_MG *)pc->data;
36549c604d5SFande Kong   PC_HMG *hmg = (PC_HMG *)mg->innerctx;
36649c604d5SFande Kong 
36749c604d5SFande Kong   PetscFunctionBegin;
36849c604d5SFande Kong   hmg->component = component;
36949c604d5SFande Kong   PetscFunctionReturn(0);
37049c604d5SFande Kong }
37149c604d5SFande Kong 
372b155ba7fSFande Kong /*@
37349c604d5SFande Kong    PCHMGSetCoarseningComponent - Set which component is used for the subspace-based coarsening algorithm
37449c604d5SFande Kong 
37549c604d5SFande Kong    Logically Collective on PC
37649c604d5SFande Kong 
37749c604d5SFande Kong    Input Parameters:
37849c604d5SFande Kong +  pc - the HMG context
37949c604d5SFande Kong -  component - which component PC will coarsen
38049c604d5SFande Kong 
38149c604d5SFande Kong    Options Database Keys:
38249c604d5SFande Kong .  -pc_hmg_coarsening_component - Which component is chosen for the subspace-based coarsening algorithm
38349c604d5SFande Kong 
38449c604d5SFande Kong    Level: beginner
38549c604d5SFande Kong 
38649c604d5SFande Kong .keywords: HMG, multigrid, interpolation, coarsening, component
38749c604d5SFande Kong 
388db781477SPatrick Sanan .seealso: `PCHMG`, `PCType`
389b155ba7fSFande Kong @*/
3909371c9d4SSatish Balay PetscErrorCode PCHMGSetCoarseningComponent(PC pc, PetscInt component) {
39149c604d5SFande Kong   PetscFunctionBegin;
39249c604d5SFande Kong   PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
393cac4c232SBarry Smith   PetscUseMethod(pc, "PCHMGSetCoarseningComponent_C", (PC, PetscInt), (pc, component));
39449c604d5SFande Kong   PetscFunctionReturn(0);
39549c604d5SFande Kong }
39649c604d5SFande Kong 
3979371c9d4SSatish Balay static PetscErrorCode PCHMGUseMatMAIJ_HMG(PC pc, PetscBool usematmaij) {
39849c604d5SFande Kong   PC_MG  *mg  = (PC_MG *)pc->data;
39949c604d5SFande Kong   PC_HMG *hmg = (PC_HMG *)mg->innerctx;
40049c604d5SFande Kong 
40149c604d5SFande Kong   PetscFunctionBegin;
40249c604d5SFande Kong   hmg->usematmaij = usematmaij;
40349c604d5SFande Kong   PetscFunctionReturn(0);
40449c604d5SFande Kong }
40549c604d5SFande Kong 
406b155ba7fSFande Kong /*@
40749c604d5SFande Kong    PCHMGUseMatMAIJ - Set a flag that indicates if or not to use MatMAIJ for interpolations for saving memory
40849c604d5SFande Kong 
40949c604d5SFande Kong    Logically Collective on PC
41049c604d5SFande Kong 
41149c604d5SFande Kong    Input Parameters:
41249c604d5SFande Kong +  pc - the HMG context
41349c604d5SFande Kong -  usematmaij - if or not to use MatMAIJ for interpolations. By default, it is true for saving memory
41449c604d5SFande Kong 
41549c604d5SFande Kong    Options Database Keys:
41649c604d5SFande Kong .  -pc_hmg_use_matmaij - <true | false >
41749c604d5SFande Kong 
41849c604d5SFande Kong    Level: beginner
41949c604d5SFande Kong 
42049c604d5SFande Kong .keywords: HMG, multigrid, interpolation, coarsening, MatMAIJ
42149c604d5SFande Kong 
422db781477SPatrick Sanan .seealso: `PCHMG`, `PCType`
423b155ba7fSFande Kong @*/
4249371c9d4SSatish Balay PetscErrorCode PCHMGUseMatMAIJ(PC pc, PetscBool usematmaij) {
42549c604d5SFande Kong   PetscFunctionBegin;
42649c604d5SFande Kong   PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
427cac4c232SBarry Smith   PetscUseMethod(pc, "PCHMGUseMatMAIJ_C", (PC, PetscBool), (pc, usematmaij));
42849c604d5SFande Kong   PetscFunctionReturn(0);
42949c604d5SFande Kong }
43049c604d5SFande Kong 
4318a2c336bSFande Kong /*MC
432fd2dd295SFande Kong    PCHMG - Hybrid of PETSc preconditioners (such as ASM, BJacobi, SOR, etc.) and Hypre BoomerAMG, GAMG or other multilevel methods. BoomerAMG, GAMG
433fd2dd295SFande Kong            or other multilevel methods is used to coarsen matrix and generate a sequence of coarse matrices and interpolations. The matrices and
434fd2dd295SFande Kong            interpolations are employed to construct PCMG, and then any available PETSc preconditioners can be chosen as smoothers and the coarse solver.
435360ee056SFande Kong 
436360ee056SFande Kong    Options Database Keys:
4378a2c336bSFande Kong +  -pc_hmg_reuse_interpolation <true | false> - Whether or not to reuse the interpolations. If true, it potentially save the compute time.
4388a2c336bSFande Kong .  -pc_hmg_use_subspace_coarsening  <true | false> - Whether or not to use subspace coarsening (that is, coarsen a submatrix).
439fd2dd295SFande Kong .  -hmg_inner_pc_type <hypre, gamg, ...> - What method is used to coarsen matrix
4404bb91820SFande Kong -  -pc_hmg_use_matmaij <true | false> - Whether or not to use MatMAIJ for multicomponent problems for saving memory
441360ee056SFande Kong 
442360ee056SFande Kong    Notes:
443360ee056SFande Kong     For multicomponent problems, we can just coarsen one submatrix associated with one particular component. In this way, the preconditioner setup
444360ee056SFande Kong     time is significantly reduced. One typical use case is neutron transport equations. There are many variables on each mesh vertex due to the
445360ee056SFande Kong     of angle and energy. Each variable, in fact, corresponds to the same PDEs but with different material properties.
446360ee056SFande Kong 
447360ee056SFande Kong    Level: beginner
448360ee056SFande Kong 
449360ee056SFande Kong     References:
450606c0280SSatish Balay .   * - Fande Kong, Yaqi Wang, Derek R Gaston, Cody J Permann, Andrew E Slaughter, Alexander D Lindsay, Richard C Martineau, A highly parallel multilevel
451360ee056SFande Kong     Newton-Krylov-Schwarz method with subspace-based coarsening and partition-based balancing for the multigroup neutron transport equations on
452360ee056SFande Kong     3D unstructured meshes, arXiv preprint arXiv:1903.03659, 2019
453360ee056SFande Kong 
454db781477SPatrick Sanan .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`, `PCMG`, `PCHYPRE`, `PCHMG`, `PCGetCoarseOperators()`, `PCGetInterpolations()`, `PCHMGSetReuseInterpolation()`, `PCHMGSetUseSubspaceCoarsening()`,
455db781477SPatrick Sanan           `PCHMGSetInnerPCType()`
456360ee056SFande Kong 
457360ee056SFande Kong M*/
4589371c9d4SSatish Balay PETSC_EXTERN PetscErrorCode PCCreate_HMG(PC pc) {
459360ee056SFande Kong   PC_HMG *hmg;
460360ee056SFande Kong   PC_MG  *mg;
461360ee056SFande Kong 
462360ee056SFande Kong   PetscFunctionBegin;
463360ee056SFande Kong   /* if type was previously mg; must manually destroy it because call to PCSetType(pc,PCMG) will not destroy it */
464dbbe0bcdSBarry Smith   PetscTryTypeMethod(pc, destroy);
4650a545947SLisandro Dalcin   pc->data = NULL;
4669566063dSJacob Faibussowitsch   PetscCall(PetscFree(((PetscObject)pc)->type_name));
467360ee056SFande Kong 
4689566063dSJacob Faibussowitsch   PetscCall(PCSetType(pc, PCMG));
4699566063dSJacob Faibussowitsch   PetscCall(PetscObjectChangeTypeName((PetscObject)pc, PCHMG));
4709566063dSJacob Faibussowitsch   PetscCall(PetscNew(&hmg));
471360ee056SFande Kong 
472360ee056SFande Kong   mg                 = (PC_MG *)pc->data;
473360ee056SFande Kong   mg->innerctx       = hmg;
474360ee056SFande Kong   hmg->reuseinterp   = PETSC_FALSE;
4758a2c336bSFande Kong   hmg->subcoarsening = PETSC_FALSE;
4764bb91820SFande Kong   hmg->usematmaij    = PETSC_TRUE;
47749c604d5SFande Kong   hmg->component     = 0;
4788a2c336bSFande Kong   hmg->innerpc       = NULL;
479360ee056SFande Kong 
480360ee056SFande Kong   pc->ops->setfromoptions = PCSetFromOptions_HMG;
481360ee056SFande Kong   pc->ops->view           = PCView_HMG;
482360ee056SFande Kong   pc->ops->destroy        = PCDestroy_HMG;
483360ee056SFande Kong   pc->ops->setup          = PCSetUp_HMG;
484fd2dd295SFande Kong 
4859566063dSJacob Faibussowitsch   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHMGSetReuseInterpolation_C", PCHMGSetReuseInterpolation_HMG));
4869566063dSJacob Faibussowitsch   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHMGSetUseSubspaceCoarsening_C", PCHMGSetUseSubspaceCoarsening_HMG));
4879566063dSJacob Faibussowitsch   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHMGSetInnerPCType_C", PCHMGSetInnerPCType_HMG));
4889566063dSJacob Faibussowitsch   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHMGSetCoarseningComponent_C", PCHMGSetCoarseningComponent_HMG));
4899566063dSJacob Faibussowitsch   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCHMGUseMatMAIJ_C", PCHMGUseMatMAIJ_HMG));
490360ee056SFande Kong   PetscFunctionReturn(0);
491360ee056SFande Kong }
492