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