Actual source code: ex31.c
1: static char help[] = "Test that shared points on interface of partitions can be rebalanced.\n\n";
2: static char FILENAME[] = "ex31.c";
4: #include <petscdmplex.h>
5: #include <petscviewerhdf5.h>
6: #include <petscsf.h>
8: typedef struct {
9: PetscBool parallel; /* Use ParMetis or Metis */
10: PetscBool useInitialGuess; /* Only active when in parallel, uses RefineKway of ParMetis */
11: PetscInt entityDepth; /* depth of the entities to rebalance ( 0 => vertices) */
12: } AppCtx;
14: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
15: {
16: options->parallel = PETSC_FALSE;
17: options->useInitialGuess = PETSC_FALSE;
18: options->entityDepth = 0;
20: PetscOptionsBegin(comm, "", "Meshing Interpolation Test Options", "DMPLEX");
21: PetscOptionsBoundedInt("-entity_depth", "Depth of the entities to rebalance (0 => vertices)", FILENAME, options->entityDepth, &options->entityDepth, NULL, 0);
22: PetscOptionsBool("-parallel", "Use ParMetis instead of Metis", FILENAME, options->parallel, &options->parallel, NULL);
23: PetscOptionsBool("-use_initial_guess", "Use RefineKway function of ParMetis", FILENAME, options->useInitialGuess, &options->useInitialGuess, NULL);
24: PetscOptionsEnd();
25: return 0;
26: }
28: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
29: {
30: DMCreate(comm, dm);
31: DMSetType(*dm, DMPLEX);
32: DMSetFromOptions(*dm);
33: DMViewFromOptions(*dm, NULL, "-dm_view");
34: return 0;
35: }
37: int main(int argc, char **argv)
38: {
39: MPI_Comm comm;
40: DM dm, dmdist;
41: PetscPartitioner part;
42: AppCtx user;
43: IS is = NULL;
44: PetscSection s = NULL, gsection = NULL;
45: PetscMPIInt size;
46: PetscSF sf;
47: PetscInt pStart, pEnd, p, minBefore, maxBefore, minAfter, maxAfter, gSizeBefore, gSizeAfter;
48: PetscBool success;
51: PetscInitialize(&argc, &argv, NULL, help);
52: comm = PETSC_COMM_WORLD;
53: MPI_Comm_size(comm, &size);
54: ProcessOptions(comm, &user);
55: CreateMesh(comm, &user, &dm);
57: /* partition dm using PETSCPARTITIONERPARMETIS */
58: DMPlexGetPartitioner(dm, &part);
59: PetscObjectSetOptionsPrefix((PetscObject)part, "p_");
60: PetscPartitionerSetType(part, PETSCPARTITIONERPARMETIS);
61: PetscPartitionerSetFromOptions(part);
62: PetscSectionCreate(comm, &s);
63: PetscPartitionerDMPlexPartition(part, dm, NULL, s, &is);
65: DMPlexDistribute(dm, 0, NULL, &dmdist);
66: if (dmdist) {
67: DMDestroy(&dm);
68: dm = dmdist;
69: }
71: /* cleanup */
72: PetscSectionDestroy(&s);
73: ISDestroy(&is);
75: /* We make a PetscSection with a DOF on every mesh entity of depth
76: * user.entityDepth, then make a global section and look at its storage size.
77: * We do the same thing after the rebalancing and then assert that the size
78: * remains the same. We also make sure that the balance has improved at least
79: * a little bit compared to the initial decomposition. */
81: if (size > 1) {
82: PetscSectionCreate(comm, &s);
83: PetscSectionSetNumFields(s, 1);
84: PetscSectionSetFieldComponents(s, 0, 1);
85: DMPlexGetDepthStratum(dm, user.entityDepth, &pStart, &pEnd);
86: PetscSectionSetChart(s, pStart, pEnd);
87: for (p = pStart; p < pEnd; ++p) {
88: PetscSectionSetDof(s, p, 1);
89: PetscSectionSetFieldDof(s, p, 0, 1);
90: }
91: PetscSectionSetUp(s);
92: DMGetPointSF(dm, &sf);
93: PetscSectionCreateGlobalSection(s, sf, PETSC_FALSE, PETSC_FALSE, &gsection);
94: PetscSectionGetStorageSize(gsection, &gSizeBefore);
95: minBefore = gSizeBefore;
96: maxBefore = gSizeBefore;
97: MPI_Allreduce(MPI_IN_PLACE, &gSizeBefore, 1, MPIU_INT, MPI_SUM, comm);
98: MPI_Allreduce(MPI_IN_PLACE, &minBefore, 1, MPIU_INT, MPI_MIN, comm);
99: MPI_Allreduce(MPI_IN_PLACE, &maxBefore, 1, MPIU_INT, MPI_MAX, comm);
100: PetscSectionDestroy(&gsection);
101: }
103: DMPlexRebalanceSharedPoints(dm, user.entityDepth, user.useInitialGuess, user.parallel, &success);
105: if (size > 1) {
106: PetscSectionCreateGlobalSection(s, sf, PETSC_FALSE, PETSC_FALSE, &gsection);
107: PetscSectionGetStorageSize(gsection, &gSizeAfter);
108: minAfter = gSizeAfter;
109: maxAfter = gSizeAfter;
110: MPI_Allreduce(MPI_IN_PLACE, &gSizeAfter, 1, MPIU_INT, MPI_SUM, comm);
111: MPI_Allreduce(MPI_IN_PLACE, &minAfter, 1, MPIU_INT, MPI_MIN, comm);
112: MPI_Allreduce(MPI_IN_PLACE, &maxAfter, 1, MPIU_INT, MPI_MAX, comm);
115: PetscSectionDestroy(&gsection);
116: PetscSectionDestroy(&s);
117: }
119: DMDestroy(&dm);
120: PetscFinalize();
121: return 0;
122: }
124: /*TEST
126: testset:
127: args: -dm_plex_dim 3 -dm_plex_simplex 0
129: test:
130: # rebalance a mesh
131: suffix: 0
132: nsize: {{2 3 4}}
133: requires: parmetis
134: args: -dm_plex_box_faces {{2,3,4 5,4,3 7,11,5}} -entity_depth {{0 1}} -parallel {{FALSE TRUE}} -use_initial_guess FALSE
136: test:
137: # rebalance a mesh but use the initial guess (uses a random algorithm and gives different results on different machines, so just check that it runs).
138: suffix: 1
139: nsize: {{2 3 4}}
140: requires: parmetis
141: args: -dm_plex_box_faces {{2,3,4 5,4,3 7,11,5}} -entity_depth {{0 1}} -parallel TRUE -use_initial_guess TRUE
143: test:
144: # no-op in serial
145: suffix: 2
146: nsize: {{1}}
147: requires: parmetis
148: args: -dm_plex_box_faces 2,3,4 -entity_depth 0 -parallel FALSE -use_initial_guess FALSE
150: TEST*/