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