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