1 static char help[] = "Performance tests for DMPlex query operations\n\n";
2
3 #include <petscdmplex.h>
4
5 typedef struct {
6 PetscInt dim; /* The topological mesh dimension */
7 PetscBool cellSimplex; /* Flag for simplices */
8 PetscBool spectral; /* Flag for spectral element layout */
9 PetscBool interpolate; /* Flag for mesh interpolation */
10 PetscReal refinementLimit; /* Maximum volume of a refined cell */
11 PetscInt numFields; /* The number of section fields */
12 PetscInt *numComponents; /* The number of field components */
13 PetscInt *numDof; /* The dof signature for the section */
14 PetscBool reuseArray; /* Pass in user allocated array to VecGetClosure() */
15 /* Test data */
16 PetscBool errors; /* Treat failures as errors */
17 PetscInt iterations; /* The number of iterations for a query */
18 PetscReal maxConeTime; /* Max time per run for DMPlexGetCone() */
19 PetscReal maxClosureTime; /* Max time per run for DMPlexGetTransitiveClosure() */
20 PetscReal maxVecClosureTime; /* Max time per run for DMPlexVecGetClosure() */
21 PetscBool printTimes; /* Print total times, do not check limits */
22 } AppCtx;
23
ProcessOptions(AppCtx * options)24 static PetscErrorCode ProcessOptions(AppCtx *options)
25 {
26 PetscInt len;
27 PetscBool flg;
28
29 PetscFunctionBegin;
30 options->dim = 2;
31 options->cellSimplex = PETSC_TRUE;
32 options->spectral = PETSC_FALSE;
33 options->interpolate = PETSC_FALSE;
34 options->refinementLimit = 0.0;
35 options->numFields = 0;
36 options->numComponents = NULL;
37 options->numDof = NULL;
38 options->reuseArray = PETSC_FALSE;
39 options->errors = PETSC_FALSE;
40 options->iterations = 1;
41 options->maxConeTime = 0.0;
42 options->maxClosureTime = 0.0;
43 options->maxVecClosureTime = 0.0;
44 options->printTimes = PETSC_FALSE;
45
46 PetscOptionsBegin(PETSC_COMM_SELF, "", "Meshing Problem Options", "DMPLEX");
47 PetscCall(PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex9.c", options->dim, &options->dim, NULL, 1, 3));
48 PetscCall(PetscOptionsBool("-cellSimplex", "Flag for simplices", "ex9.c", options->cellSimplex, &options->cellSimplex, NULL));
49 PetscCall(PetscOptionsBool("-spectral", "Flag for spectral element layout", "ex9.c", options->spectral, &options->spectral, NULL));
50 PetscCall(PetscOptionsBool("-interpolate", "Flag for mesh interpolation", "ex9.c", options->interpolate, &options->interpolate, NULL));
51 PetscCall(PetscOptionsReal("-refinement_limit", "The maximum volume of a refined cell", "ex9.c", options->refinementLimit, &options->refinementLimit, NULL));
52 PetscCall(PetscOptionsBoundedInt("-num_fields", "The number of section fields", "ex9.c", options->numFields, &options->numFields, NULL, 0));
53 if (options->numFields) {
54 len = options->numFields;
55 PetscCall(PetscMalloc1(len, &options->numComponents));
56 PetscCall(PetscOptionsIntArray("-num_components", "The number of components per field", "ex9.c", options->numComponents, &len, &flg));
57 PetscCheck(!flg || !(len != options->numFields), PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Length of components array is %" PetscInt_FMT " should be %" PetscInt_FMT, len, options->numFields);
58 }
59 len = (options->dim + 1) * PetscMax(1, options->numFields);
60 PetscCall(PetscMalloc1(len, &options->numDof));
61 PetscCall(PetscOptionsIntArray("-num_dof", "The dof signature for the section", "ex9.c", options->numDof, &len, &flg));
62 PetscCheck(!flg || len == (options->dim + 1) * PetscMax(1, options->numFields), PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Length of dof array is %" PetscInt_FMT " should be %" PetscInt_FMT, len, (options->dim + 1) * PetscMax(1, options->numFields));
63
64 /* We are specifying the scalar dof, so augment it for multiple components */
65 {
66 PetscInt f, d;
67
68 for (f = 0; f < options->numFields; ++f) {
69 for (d = 0; d <= options->dim; ++d) options->numDof[f * (options->dim + 1) + d] *= options->numComponents[f];
70 }
71 }
72
73 PetscCall(PetscOptionsBool("-reuse_array", "Pass in user allocated array to VecGetClosure()", "ex9.c", options->reuseArray, &options->reuseArray, NULL));
74 PetscCall(PetscOptionsBool("-errors", "Treat failures as errors", "ex9.c", options->errors, &options->errors, NULL));
75 PetscCall(PetscOptionsBoundedInt("-iterations", "The number of iterations for a query", "ex9.c", options->iterations, &options->iterations, NULL, 0));
76 PetscCall(PetscOptionsReal("-max_cone_time", "The maximum time per run for DMPlexGetCone()", "ex9.c", options->maxConeTime, &options->maxConeTime, NULL));
77 PetscCall(PetscOptionsReal("-max_closure_time", "The maximum time per run for DMPlexGetTransitiveClosure()", "ex9.c", options->maxClosureTime, &options->maxClosureTime, NULL));
78 PetscCall(PetscOptionsReal("-max_vec_closure_time", "The maximum time per run for DMPlexVecGetClosure()", "ex9.c", options->maxVecClosureTime, &options->maxVecClosureTime, NULL));
79 PetscCall(PetscOptionsBool("-print_times", "Print total times, do not check limits", "ex9.c", options->printTimes, &options->printTimes, NULL));
80 PetscOptionsEnd();
81 PetscFunctionReturn(PETSC_SUCCESS);
82 }
83
CreateSimplex_2D(MPI_Comm comm,DM * newdm)84 static PetscErrorCode CreateSimplex_2D(MPI_Comm comm, DM *newdm)
85 {
86 DM dm;
87 PetscInt numPoints[2] = {4, 2};
88 PetscInt coneSize[6] = {3, 3, 0, 0, 0, 0};
89 PetscInt cones[6] = {2, 3, 4, 5, 4, 3};
90 PetscInt coneOrientations[6] = {0, 0, 0, 0, 0, 0};
91 PetscScalar vertexCoords[8] = {-0.5, 0.5, 0.0, 0.0, 0.0, 1.0, 0.5, 0.5};
92 PetscInt markerPoints[8] = {2, 1, 3, 1, 4, 1, 5, 1};
93 PetscInt dim = 2, depth = 1, p;
94
95 PetscFunctionBegin;
96 PetscCall(DMCreate(comm, &dm));
97 PetscCall(PetscObjectSetName((PetscObject)dm, "triangular"));
98 PetscCall(DMSetType(dm, DMPLEX));
99 PetscCall(DMSetDimension(dm, dim));
100 PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords));
101 for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
102 *newdm = dm;
103 PetscFunctionReturn(PETSC_SUCCESS);
104 }
105
CreateSimplex_3D(MPI_Comm comm,DM * newdm)106 static PetscErrorCode CreateSimplex_3D(MPI_Comm comm, DM *newdm)
107 {
108 DM dm;
109 PetscInt numPoints[2] = {5, 2};
110 PetscInt coneSize[23] = {4, 4, 0, 0, 0, 0, 0};
111 PetscInt cones[8] = {2, 4, 3, 5, 3, 4, 6, 5};
112 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
113 PetscScalar vertexCoords[15] = {0.0, 0.0, -0.5, 0.0, -0.5, 0.0, 1.0, 0.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0, 0.5};
114 PetscInt markerPoints[10] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1};
115 PetscInt dim = 3, depth = 1, p;
116
117 PetscFunctionBegin;
118 PetscCall(DMCreate(comm, &dm));
119 PetscCall(PetscObjectSetName((PetscObject)dm, "tetrahedral"));
120 PetscCall(DMSetType(dm, DMPLEX));
121 PetscCall(DMSetDimension(dm, dim));
122 PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords));
123 for (p = 0; p < 5; ++p) PetscCall(DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
124 *newdm = dm;
125 PetscFunctionReturn(PETSC_SUCCESS);
126 }
127
CreateQuad_2D(MPI_Comm comm,DM * newdm)128 static PetscErrorCode CreateQuad_2D(MPI_Comm comm, DM *newdm)
129 {
130 DM dm;
131 PetscInt numPoints[2] = {6, 2};
132 PetscInt coneSize[8] = {4, 4, 0, 0, 0, 0, 0, 0};
133 PetscInt cones[8] = {2, 3, 4, 5, 3, 6, 7, 4};
134 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
135 PetscScalar vertexCoords[12] = {-0.5, 0.0, 0.0, 0.0, 0.0, 1.0, -0.5, 1.0, 0.5, 0.0, 0.5, 1.0};
136 PetscInt markerPoints[12] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1, 7, 1};
137 PetscInt dim = 2, depth = 1, p;
138
139 PetscFunctionBegin;
140 PetscCall(DMCreate(comm, &dm));
141 PetscCall(PetscObjectSetName((PetscObject)dm, "quadrilateral"));
142 PetscCall(DMSetType(dm, DMPLEX));
143 PetscCall(DMSetDimension(dm, dim));
144 PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords));
145 for (p = 0; p < 6; ++p) PetscCall(DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
146 *newdm = dm;
147 PetscFunctionReturn(PETSC_SUCCESS);
148 }
149
CreateHex_3D(MPI_Comm comm,DM * newdm)150 static PetscErrorCode CreateHex_3D(MPI_Comm comm, DM *newdm)
151 {
152 DM dm;
153 PetscInt numPoints[2] = {12, 2};
154 PetscInt coneSize[14] = {8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
155 PetscInt cones[16] = {2, 5, 4, 3, 6, 7, 8, 9, 3, 4, 11, 10, 7, 12, 13, 8};
156 PetscInt coneOrientations[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
157 PetscScalar vertexCoords[36] = {-0.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, -0.5, 1.0, 0.0, -0.5, 0.0, 1.0, 0.0, 0.0, 1.0, 0.0, 1.0, 1.0, -0.5, 1.0, 1.0, 0.5, 0.0, 0.0, 0.5, 1.0, 0.0, 0.5, 0.0, 1.0, 0.5, 1.0, 1.0};
158 PetscInt markerPoints[24] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1, 7, 1, 8, 1, 9, 1, 10, 1, 11, 1, 12, 1, 13, 1};
159 PetscInt dim = 3, depth = 1, p;
160
161 PetscFunctionBegin;
162 PetscCall(DMCreate(comm, &dm));
163 PetscCall(PetscObjectSetName((PetscObject)dm, "hexahedral"));
164 PetscCall(DMSetType(dm, DMPLEX));
165 PetscCall(DMSetDimension(dm, dim));
166 PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords));
167 for (p = 0; p < 12; ++p) PetscCall(DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
168 *newdm = dm;
169 PetscFunctionReturn(PETSC_SUCCESS);
170 }
171
CreateMesh(MPI_Comm comm,AppCtx * user,DM * newdm)172 static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *newdm)
173 {
174 PetscInt dim = user->dim;
175 PetscBool cellSimplex = user->cellSimplex;
176
177 PetscFunctionBegin;
178 switch (dim) {
179 case 2:
180 if (cellSimplex) {
181 PetscCall(CreateSimplex_2D(comm, newdm));
182 } else {
183 PetscCall(CreateQuad_2D(comm, newdm));
184 }
185 break;
186 case 3:
187 if (cellSimplex) {
188 PetscCall(CreateSimplex_3D(comm, newdm));
189 } else {
190 PetscCall(CreateHex_3D(comm, newdm));
191 }
192 break;
193 default:
194 SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make meshes for dimension %" PetscInt_FMT, dim);
195 }
196 if (user->refinementLimit > 0.0) {
197 DM rdm;
198 const char *name;
199
200 PetscCall(DMPlexSetRefinementUniform(*newdm, PETSC_FALSE));
201 PetscCall(DMPlexSetRefinementLimit(*newdm, user->refinementLimit));
202 PetscCall(DMRefine(*newdm, PETSC_COMM_SELF, &rdm));
203 PetscCall(PetscObjectGetName((PetscObject)*newdm, &name));
204 PetscCall(PetscObjectSetName((PetscObject)rdm, name));
205 PetscCall(DMDestroy(newdm));
206 *newdm = rdm;
207 }
208 if (user->interpolate) {
209 DM idm;
210
211 PetscCall(DMPlexInterpolate(*newdm, &idm));
212 PetscCall(DMDestroy(newdm));
213 *newdm = idm;
214 }
215 PetscCall(DMSetFromOptions(*newdm));
216 PetscFunctionReturn(PETSC_SUCCESS);
217 }
218
TestCone(DM dm,AppCtx * user)219 static PetscErrorCode TestCone(DM dm, AppCtx *user)
220 {
221 PetscInt numRuns, cStart, cEnd, c, i;
222 PetscReal maxTimePerRun = user->maxConeTime;
223 PetscLogStage stage;
224 PetscLogEvent event;
225 PetscEventPerfInfo eventInfo;
226 MPI_Comm comm;
227 PetscMPIInt rank;
228
229 PetscFunctionBegin;
230 PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
231 PetscCallMPI(MPI_Comm_rank(comm, &rank));
232 PetscCall(PetscLogStageRegister("DMPlex Cone Test", &stage));
233 PetscCall(PetscLogEventRegister("Cone", PETSC_OBJECT_CLASSID, &event));
234 PetscCall(PetscLogStagePush(stage));
235 PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd));
236 PetscCall(PetscLogEventBegin(event, 0, 0, 0, 0));
237 for (i = 0; i < user->iterations; ++i) {
238 for (c = cStart; c < cEnd; ++c) {
239 const PetscInt *cone;
240
241 PetscCall(DMPlexGetCone(dm, c, &cone));
242 }
243 }
244 PetscCall(PetscLogEventEnd(event, 0, 0, 0, 0));
245 PetscCall(PetscLogStagePop());
246
247 PetscCall(PetscLogEventGetPerfInfo(stage, event, &eventInfo));
248 numRuns = (cEnd - cStart) * user->iterations;
249 PetscCheck(eventInfo.count == 1, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Number of event calls %d should be 1", eventInfo.count);
250 PetscCheck((PetscInt)eventInfo.flops == 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Number of event flops %" PetscInt_FMT " should be 0", (PetscInt)eventInfo.flops);
251 if (user->printTimes) {
252 PetscCall(PetscSynchronizedPrintf(comm, "[%d] Cones: %" PetscInt_FMT " Total time: %.3es Average time per cone: %.3es\n", rank, numRuns, eventInfo.time, eventInfo.time / numRuns));
253 PetscCall(PetscSynchronizedFlush(comm, PETSC_STDOUT));
254 } else if (eventInfo.time > maxTimePerRun * numRuns) {
255 PetscCall(PetscSynchronizedPrintf(comm, "[%d] Cones: %" PetscInt_FMT " Average time per cone: %gs standard: %gs\n", rank, numRuns, eventInfo.time / numRuns, (double)maxTimePerRun));
256 PetscCall(PetscSynchronizedFlush(comm, PETSC_STDOUT));
257 PetscCheck(!user->errors, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Average time for cone %g > standard %g", eventInfo.time / numRuns, (double)maxTimePerRun);
258 }
259 PetscFunctionReturn(PETSC_SUCCESS);
260 }
261
TestTransitiveClosure(DM dm,AppCtx * user)262 static PetscErrorCode TestTransitiveClosure(DM dm, AppCtx *user)
263 {
264 PetscInt numRuns, cStart, cEnd, c, i;
265 PetscReal maxTimePerRun = user->maxClosureTime;
266 PetscLogStage stage;
267 PetscLogEvent event;
268 PetscEventPerfInfo eventInfo;
269 MPI_Comm comm;
270 PetscMPIInt rank;
271
272 PetscFunctionBegin;
273 PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
274 PetscCallMPI(MPI_Comm_rank(comm, &rank));
275 PetscCall(PetscLogStageRegister("DMPlex Transitive Closure Test", &stage));
276 PetscCall(PetscLogEventRegister("TransitiveClosure", PETSC_OBJECT_CLASSID, &event));
277 PetscCall(PetscLogStagePush(stage));
278 PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd));
279 PetscCall(PetscLogEventBegin(event, 0, 0, 0, 0));
280 for (i = 0; i < user->iterations; ++i) {
281 for (c = cStart; c < cEnd; ++c) {
282 PetscInt *closure = NULL;
283 PetscInt closureSize;
284
285 PetscCall(DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure));
286 PetscCall(DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure));
287 }
288 }
289 PetscCall(PetscLogEventEnd(event, 0, 0, 0, 0));
290 PetscCall(PetscLogStagePop());
291
292 PetscCall(PetscLogEventGetPerfInfo(stage, event, &eventInfo));
293 numRuns = (cEnd - cStart) * user->iterations;
294 PetscCheck(eventInfo.count == 1, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Number of event calls %d should be 1", eventInfo.count);
295 PetscCheck((PetscInt)eventInfo.flops == 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Number of event flops %" PetscInt_FMT " should be 0", (PetscInt)eventInfo.flops);
296 if (user->printTimes) {
297 PetscCall(PetscSynchronizedPrintf(comm, "[%d] Closures: %" PetscInt_FMT " Total time: %.3es Average time per cone: %.3es\n", rank, numRuns, eventInfo.time, eventInfo.time / numRuns));
298 PetscCall(PetscSynchronizedFlush(comm, PETSC_STDOUT));
299 } else if (eventInfo.time > maxTimePerRun * numRuns) {
300 PetscCall(PetscSynchronizedPrintf(comm, "[%d] Closures: %" PetscInt_FMT " Average time per cone: %gs standard: %gs\n", rank, numRuns, eventInfo.time / numRuns, (double)maxTimePerRun));
301 PetscCall(PetscSynchronizedFlush(comm, PETSC_STDOUT));
302 PetscCheck(!user->errors, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Average time for closure %g > standard %g", eventInfo.time / numRuns, (double)maxTimePerRun);
303 }
304 PetscFunctionReturn(PETSC_SUCCESS);
305 }
306
TestVecClosure(DM dm,PetscBool useIndex,PetscBool useSpectral,AppCtx * user)307 static PetscErrorCode TestVecClosure(DM dm, PetscBool useIndex, PetscBool useSpectral, AppCtx *user)
308 {
309 PetscSection s;
310 Vec v;
311 PetscInt numRuns, cStart, cEnd, c, i;
312 PetscScalar tmpArray[64];
313 PetscScalar *userArray = user->reuseArray ? tmpArray : NULL;
314 PetscReal maxTimePerRun = user->maxVecClosureTime;
315 PetscLogStage stage;
316 PetscLogEvent event;
317 PetscEventPerfInfo eventInfo;
318 MPI_Comm comm;
319 PetscMPIInt rank;
320
321 PetscFunctionBegin;
322 PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
323 PetscCallMPI(MPI_Comm_rank(comm, &rank));
324 if (useIndex) {
325 if (useSpectral) {
326 PetscCall(PetscLogStageRegister("DMPlex Vector Closure with Index Test", &stage));
327 PetscCall(PetscLogEventRegister("VecClosureInd", PETSC_OBJECT_CLASSID, &event));
328 } else {
329 PetscCall(PetscLogStageRegister("DMPlex Vector Spectral Closure with Index Test", &stage));
330 PetscCall(PetscLogEventRegister("VecClosureSpecInd", PETSC_OBJECT_CLASSID, &event));
331 }
332 } else {
333 if (useSpectral) {
334 PetscCall(PetscLogStageRegister("DMPlex Vector Spectral Closure Test", &stage));
335 PetscCall(PetscLogEventRegister("VecClosureSpec", PETSC_OBJECT_CLASSID, &event));
336 } else {
337 PetscCall(PetscLogStageRegister("DMPlex Vector Closure Test", &stage));
338 PetscCall(PetscLogEventRegister("VecClosure", PETSC_OBJECT_CLASSID, &event));
339 }
340 }
341 PetscCall(PetscLogStagePush(stage));
342 PetscCall(DMSetNumFields(dm, user->numFields));
343 PetscCall(DMPlexCreateSection(dm, NULL, user->numComponents, user->numDof, 0, NULL, NULL, NULL, NULL, &s));
344 PetscCall(DMSetLocalSection(dm, s));
345 if (useIndex) PetscCall(DMPlexCreateClosureIndex(dm, s));
346 if (useSpectral) PetscCall(DMPlexSetClosurePermutationTensor(dm, PETSC_DETERMINE, s));
347 PetscCall(PetscSectionDestroy(&s));
348 PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd));
349 PetscCall(DMGetLocalVector(dm, &v));
350 PetscCall(PetscLogEventBegin(event, 0, 0, 0, 0));
351 for (i = 0; i < user->iterations; ++i) {
352 for (c = cStart; c < cEnd; ++c) {
353 PetscScalar *closure = userArray;
354 PetscInt closureSize = 64;
355
356 PetscCall(DMPlexVecGetClosure(dm, s, v, c, &closureSize, &closure));
357 if (!user->reuseArray) PetscCall(DMPlexVecRestoreClosure(dm, s, v, c, &closureSize, &closure));
358 }
359 }
360 PetscCall(PetscLogEventEnd(event, 0, 0, 0, 0));
361 PetscCall(DMRestoreLocalVector(dm, &v));
362 PetscCall(PetscLogStagePop());
363
364 PetscCall(PetscLogEventGetPerfInfo(stage, event, &eventInfo));
365 numRuns = (cEnd - cStart) * user->iterations;
366 PetscCheck(eventInfo.count == 1, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Number of event calls %d should be 1", eventInfo.count);
367 PetscCheck((PetscInt)eventInfo.flops == 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Number of event flops %" PetscInt_FMT " should be 0", (PetscInt)eventInfo.flops);
368 if (user->printTimes || eventInfo.time > maxTimePerRun * numRuns) {
369 const char *title = "VecClosures";
370 const char *titleIndex = "VecClosures with Index";
371 const char *titleSpec = "VecClosures Spectral";
372 const char *titleSpecIndex = "VecClosures Spectral with Index";
373
374 if (user->printTimes) {
375 PetscCall(PetscSynchronizedPrintf(comm, "[%d] %s: %" PetscInt_FMT " Total time: %.3es Average time per vector closure: %.3es\n", rank, useIndex ? (useSpectral ? titleSpecIndex : titleIndex) : (useSpectral ? titleSpec : title), numRuns,
376 eventInfo.time, eventInfo.time / numRuns));
377 PetscCall(PetscSynchronizedFlush(comm, PETSC_STDOUT));
378 } else {
379 PetscCall(PetscSynchronizedPrintf(comm, "[%d] %s: %" PetscInt_FMT " Average time per vector closure: %gs standard: %gs\n", rank, useIndex ? (useSpectral ? titleSpecIndex : titleIndex) : (useSpectral ? titleSpec : title), numRuns, eventInfo.time / numRuns, (double)maxTimePerRun));
380 PetscCall(PetscSynchronizedFlush(comm, PETSC_STDOUT));
381 PetscCheck(!user->errors, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Average time for vector closure %g > standard %g", eventInfo.time / numRuns, (double)maxTimePerRun);
382 }
383 }
384 PetscFunctionReturn(PETSC_SUCCESS);
385 }
386
CleanupContext(AppCtx * user)387 static PetscErrorCode CleanupContext(AppCtx *user)
388 {
389 PetscFunctionBegin;
390 PetscCall(PetscFree(user->numComponents));
391 PetscCall(PetscFree(user->numDof));
392 PetscFunctionReturn(PETSC_SUCCESS);
393 }
394
main(int argc,char ** argv)395 int main(int argc, char **argv)
396 {
397 DM dm;
398 AppCtx user;
399
400 PetscFunctionBeginUser;
401 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
402 PetscCall(ProcessOptions(&user));
403 PetscCall(PetscLogDefaultBegin());
404 PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm));
405 PetscCall(TestCone(dm, &user));
406 PetscCall(TestTransitiveClosure(dm, &user));
407 PetscCall(TestVecClosure(dm, PETSC_FALSE, PETSC_FALSE, &user));
408 PetscCall(TestVecClosure(dm, PETSC_TRUE, PETSC_FALSE, &user));
409 if (!user.cellSimplex && user.spectral) {
410 PetscCall(TestVecClosure(dm, PETSC_FALSE, PETSC_TRUE, &user));
411 PetscCall(TestVecClosure(dm, PETSC_TRUE, PETSC_TRUE, &user));
412 }
413 PetscCall(DMDestroy(&dm));
414 PetscCall(CleanupContext(&user));
415 PetscCall(PetscFinalize());
416 return 0;
417 }
418
419 /*TEST
420
421 build:
422 requires: defined(PETSC_USE_LOG)
423
424 # 2D Simplex P_1 scalar tests
425 testset:
426 args: -num_dof 1,0,0 -iterations 2 -print_times
427 test:
428 suffix: correctness_0
429 test:
430 suffix: correctness_1
431 args: -interpolate -dm_refine 2
432 test:
433 suffix: correctness_2
434 requires: triangle
435 args: -interpolate -dm_refine 5
436 test:
437 suffix: 0
438 TODO: Only for performance testing
439 args: -num_dof 1,0,0 -iterations 10000 -max_cone_time 1.1e-8 -max_closure_time 1.3e-7 -max_vec_closure_time 3.6e-7
440 test:
441 suffix: 1
442 requires: triangle
443 TODO: Only for performance testing
444 args: -refinement_limit 1.0e-5 -num_dof 1,0,0 -iterations 2 -max_cone_time 2.1e-8 -max_closure_time 1.5e-7 -max_vec_closure_time 3.6e-7
445 test:
446 suffix: 2
447 TODO: Only for performance testing
448 args: -num_fields 1 -num_components 1 -num_dof 1,0,0 -iterations 10000 -max_cone_time 1.1e-8 -max_closure_time 1.3e-7 -max_vec_closure_time 4.5e-7
449 test:
450 suffix: 3
451 requires: triangle
452 TODO: Only for performance testing
453 args: -refinement_limit 1.0e-5 -num_fields 1 -num_components 1 -num_dof 1,0,0 -iterations 2 -max_cone_time 2.1e-8 -max_closure_time 1.5e-7 -max_vec_closure_time 4.7e-7
454 test:
455 suffix: 4
456 TODO: Only for performance testing
457 args: -interpolate -num_dof 1,0,0 -iterations 10000 -max_cone_time 1.1e-8 -max_closure_time 6.5e-7 -max_vec_closure_time 1.0e-6
458 test:
459 suffix: 5
460 requires: triangle
461 TODO: Only for performance testing
462 args: -interpolate -refinement_limit 1.0e-4 -num_dof 1,0,0 -iterations 2 -max_cone_time 2.1e-8 -max_closure_time 6.5e-7 -max_vec_closure_time 1.0e-6
463 test:
464 suffix: 6
465 TODO: Only for performance testing
466 args: -interpolate -num_fields 1 -num_components 1 -num_dof 1,0,0 -iterations 10000 -max_cone_time 1.1e-8 -max_closure_time 6.5e-7 -max_vec_closure_time 1.1e-6
467 test:
468 suffix: 7
469 requires: triangle
470 TODO: Only for performance testing
471 args: -interpolate -refinement_limit 1.0e-4 -num_fields 1 -num_components 1 -num_dof 1,0,0 -iterations 2 -max_cone_time 2.1e-8 -max_closure_time 6.5e-7 -max_vec_closure_time 1.2e-6
472
473 # 2D Simplex P_1 vector tests
474 # 2D Simplex P_2 scalar tests
475 # 2D Simplex P_2 vector tests
476 # 2D Simplex P_2/P_1 vector/scalar tests
477 # 2D Quad P_1 scalar tests
478 # 2D Quad P_1 vector tests
479 # 2D Quad P_2 scalar tests
480 # 2D Quad P_2 vector tests
481 # 3D Simplex P_1 scalar tests
482 # 3D Simplex P_1 vector tests
483 # 3D Simplex P_2 scalar tests
484 # 3D Simplex P_2 vector tests
485 # 3D Hex P_1 scalar tests
486 # 3D Hex P_1 vector tests
487 # 3D Hex P_2 scalar tests
488 # 3D Hex P_2 vector tests
489
490 TEST*/
491