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