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