Actual source code: ex10.c

  1: #include <petscksp.h>
  2: #include <petscpc.h>
  3: #include <petscviewer.h>

  5: typedef struct {
  6:   PetscInt  num_levels;
  7:   PetscInt *n_per_level;
  8:   Mat       stiff;
  9:   Mat      *ProlongationOps;
 10:   PetscBT  *CFMarkers;
 11:   KSP       kspHypre;
 12: } *DataCompression;

 14: PetscErrorCode Create1dLaplacian(PetscInt, Mat *);
 15: PetscErrorCode DataCompExportMats(DataCompression);
 16: PetscErrorCode DataCompDestroy(DataCompression);

 18: int main(int Argc, char **Args)
 19: {
 20:   PetscInt        n_nodes = 33;
 21:   Vec             x, b;
 22:   PC              pcHypre;
 23:   DataCompression data_comp;

 25:   PetscFunctionBeginUser;
 26:   PetscCall(PetscInitialize(&Argc, &Args, NULL, NULL));

 28:   PetscCall(PetscNew(&data_comp));

 30:   // Creating stiffness matrix
 31:   PetscCall(Create1dLaplacian(n_nodes, &data_comp->stiff));
 32:   PetscCall(PetscObjectSetName((PetscObject)data_comp->stiff, "Stiffness"));

 34:   // Set-up BoomerAMG PC to get Prolongation Operators and Coarse/Fine splittings
 35:   PetscCall(KSPCreate(PETSC_COMM_WORLD, &data_comp->kspHypre));
 36:   PetscCall(KSPSetType(data_comp->kspHypre, KSPRICHARDSON));
 37:   PetscCall(KSPGetPC(data_comp->kspHypre, &pcHypre));
 38:   PetscCall(PCSetType(pcHypre, PCHYPRE));
 39:   PetscCall(PCHYPRESetType(pcHypre, "boomeramg"));
 40:   PetscCall(PCSetFromOptions(pcHypre));
 41:   PetscCall(PCSetOperators(pcHypre, data_comp->stiff, data_comp->stiff));
 42:   PetscCall(PCSetUp(pcHypre));

 44:   PetscCall(MatCreateVecs(data_comp->stiff, &x, &b));
 45:   PetscCall(PCApply(pcHypre, x, b));
 46:   PetscCall(VecDestroy(&x));
 47:   PetscCall(VecDestroy(&b));

 49:   //Viewing the PC and Extracting the Prolongation Operators and CFMarkers
 50:   PetscCall(PCView(pcHypre, NULL));
 51:   PetscCall(PCGetInterpolations(pcHypre, &data_comp->num_levels, &data_comp->ProlongationOps));
 52:   PetscCall(PCHYPREGetCFMarkers(pcHypre, &data_comp->n_per_level, &data_comp->CFMarkers));

 54:   PetscCall(DataCompExportMats(data_comp));
 55:   PetscCall(DataCompDestroy(data_comp));
 56:   PetscCall(PetscFinalize());
 57:   return 0;
 58: }

 60: PetscErrorCode Create1dLaplacian(PetscInt n, Mat *mat)
 61: {
 62:   PetscFunctionBeginUser;
 63:   PetscCall(MatCreateSeqAIJ(PETSC_COMM_SELF, n, n, 3, NULL, mat));
 64:   PetscCall(MatSetValue(*mat, n - 1, n - 1, 2.0, INSERT_VALUES));
 65:   for (PetscInt i = 0; i < n - 1; i++) {
 66:     PetscCall(MatSetValue(*mat, i, i, 2.0, INSERT_VALUES));
 67:     PetscCall(MatSetValue(*mat, i + 1, i, -1.0, INSERT_VALUES));
 68:     PetscCall(MatSetValue(*mat, i, i + 1, -1.0, INSERT_VALUES));
 69:   }
 70:   PetscCall(MatAssemblyBegin(*mat, MAT_FINAL_ASSEMBLY));
 71:   PetscCall(MatAssemblyEnd(*mat, MAT_FINAL_ASSEMBLY));
 72:   PetscFunctionReturn(PETSC_SUCCESS);
 73: }

 75: PetscErrorCode DataCompExportMats(DataCompression data_comp)
 76: {
 77:   PetscFunctionBeginUser;
 78:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Num levels: %" PetscInt_FMT "\n", data_comp->num_levels));
 79:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, " -- Nodes per level --\n"));
 80:   for (PetscInt i = 0; i < data_comp->num_levels; i++) { PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Level %" PetscInt_FMT ": %" PetscInt_FMT "\n", i, data_comp->n_per_level[i])); }

 82:   for (PetscInt i = 0; i < data_comp->num_levels - 1; i++) {
 83:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Prolongation Operator - Level %" PetscInt_FMT "\n", i));
 84:     PetscCall(PetscObjectSetName((PetscObject)data_comp->ProlongationOps[i], "P"));
 85:     PetscCall(MatView(data_comp->ProlongationOps[i], PETSC_VIEWER_STDOUT_WORLD));
 86:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\n"));
 87:   }

 89:   for (PetscInt i = 0; i < data_comp->num_levels - 1; i++) {
 90:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Coarse/Fine splitting - Level %" PetscInt_FMT "\n", i + 1));
 91:     PetscCall(PetscBTView(data_comp->n_per_level[i + 1], data_comp->CFMarkers[i], PETSC_VIEWER_STDOUT_WORLD));
 92:   }
 93:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Stiffness matrix, sparse format:\n"));
 94:   PetscCall(MatViewFromOptions(data_comp->stiff, NULL, "-mat_view_stiff"));
 95:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Finished calling the Viewer functions\n"));
 96:   PetscFunctionReturn(PETSC_SUCCESS);
 97: }

 99: PetscErrorCode DataCompDestroy(DataCompression data_comp)
100: {
101:   PetscFunctionBeginUser;
102:   if (data_comp == NULL) PetscFunctionReturn(PETSC_SUCCESS);
103:   PetscCall(MatDestroy(&data_comp->stiff));
104:   PetscCall(KSPDestroy(&data_comp->kspHypre));
105:   for (PetscInt i = 0; i < data_comp->num_levels - 1; i++) {
106:     PetscCall(MatDestroy(&data_comp->ProlongationOps[i]));
107:     PetscCall(PetscBTDestroy(&data_comp->CFMarkers[i]));
108:   }
109:   PetscCall(PetscFree(data_comp->ProlongationOps));
110:   PetscCall(PetscFree(data_comp->n_per_level));
111:   PetscCall(PetscFree(data_comp->CFMarkers));
112:   PetscCall(PetscFree(data_comp));
113:   PetscFunctionReturn(PETSC_SUCCESS);
114: }

116: /*TEST

118:    test:
119:       requires: hypre !complex !defined(PETSC_HAVE_HYPRE_DEVICE)
120:       args: -pc_hypre_boomeramg_coarsen_type modifiedRuge-Stueben -pc_hypre_boomeramg_interp_type classical -pc_hypre_boomeramg_strong_threshold 0.25 pc_hypre_boomeramg_numfunctions 1 -pc_hypre_boomeramg_max_row_sum 1.0 -mat_view_stiff

122: TEST*/