1 2 #include <petscfe.h> 3 #include <petscdmplex.h> 4 #include <petsc/private/hashmap.h> 5 #include <petsc/private/dmpleximpl.h> 6 #include <petsc/private/petscfeimpl.h> 7 8 const char help[] = "Test PETSCDUALSPACELAGRANGE\n"; 9 10 typedef struct _PetscHashLagKey { 11 PetscInt dim; 12 PetscInt order; 13 PetscInt formDegree; 14 PetscBool trimmed; 15 PetscInt tensor; 16 PetscBool continuous; 17 } PetscHashLagKey; 18 19 #define PetscHashLagKeyHash(key) \ 20 PetscHashCombine(PetscHashCombine(PetscHashCombine(PetscHashInt((key).dim), PetscHashInt((key).order)), PetscHashInt((key).formDegree)), PetscHashCombine(PetscHashCombine(PetscHashInt((key).trimmed), PetscHashInt((key).tensor)), PetscHashInt((key).continuous))) 21 22 #define PetscHashLagKeyEqual(k1, k2) \ 23 (((k1).dim == (k2).dim) ? ((k1).order == (k2).order) ? ((k1).formDegree == (k2).formDegree) ? ((k1).trimmed == (k2).trimmed) ? ((k1).tensor == (k2).tensor) ? ((k1).continuous == (k2).continuous) : 0 : 0 : 0 : 0 : 0) 24 25 PETSC_HASH_MAP(HashLag, PetscHashLagKey, PetscInt, PetscHashLagKeyHash, PetscHashLagKeyEqual, 0) 26 27 static PetscErrorCode ExpectedNumDofs_Total(PetscInt dim, PetscInt order, PetscInt formDegree, PetscBool trimmed, PetscInt tensor, PetscInt nCopies, PetscInt *nDofs); 28 static PetscErrorCode ExpectedNumDofs_Interior(PetscInt dim, PetscInt order, PetscInt formDegree, PetscBool trimmed, PetscInt tensor, PetscInt nCopies, PetscInt *nDofs); 29 30 static PetscErrorCode ExpectedNumDofs_Total(PetscInt dim, PetscInt order, PetscInt formDegree, PetscBool trimmed, PetscInt tensor, PetscInt nCopies, PetscInt *nDofs) 31 { 32 PetscFunctionBegin; 33 formDegree = PetscAbsInt(formDegree); 34 /* see femtable.org for the source of most of these values */ 35 *nDofs = -1; 36 if (tensor == 0) { /* simplex spaces */ 37 if (trimmed) { 38 PetscInt rnchooserk; 39 PetscInt rkm1choosek; 40 41 PetscCall(PetscDTBinomialInt(order + dim, order + formDegree, &rnchooserk)); 42 PetscCall(PetscDTBinomialInt(order + formDegree - 1, formDegree, &rkm1choosek)); 43 *nDofs = rnchooserk * rkm1choosek * nCopies; 44 } else { 45 PetscInt rnchooserk; 46 PetscInt rkchoosek; 47 48 PetscCall(PetscDTBinomialInt(order + dim, order + formDegree, &rnchooserk)); 49 PetscCall(PetscDTBinomialInt(order + formDegree, formDegree, &rkchoosek)); 50 *nDofs = rnchooserk * rkchoosek * nCopies; 51 } 52 } else if (tensor == 1) { /* hypercubes */ 53 if (trimmed) { 54 PetscInt nchoosek; 55 PetscInt rpowk, rp1pownmk; 56 57 PetscCall(PetscDTBinomialInt(dim, formDegree, &nchoosek)); 58 rpowk = PetscPowInt(order, formDegree); 59 rp1pownmk = PetscPowInt(order + 1, dim - formDegree); 60 *nDofs = nchoosek * rpowk * rp1pownmk * nCopies; 61 } else { 62 PetscInt nchoosek; 63 PetscInt rp1pown; 64 65 PetscCall(PetscDTBinomialInt(dim, formDegree, &nchoosek)); 66 rp1pown = PetscPowInt(order + 1, dim); 67 *nDofs = nchoosek * rp1pown * nCopies; 68 } 69 } else { /* prism */ 70 PetscInt tracek = 0; 71 PetscInt tracekm1 = 0; 72 PetscInt fiber0 = 0; 73 PetscInt fiber1 = 0; 74 75 if (formDegree < dim) { 76 PetscCall(ExpectedNumDofs_Total(dim - 1, order, formDegree, trimmed, 0, 1, &tracek)); 77 PetscCall(ExpectedNumDofs_Total(1, order, 0, trimmed, 0, 1, &fiber0)); 78 } 79 if (formDegree > 0) { 80 PetscCall(ExpectedNumDofs_Total(dim - 1, order, formDegree - 1, trimmed, 0, 1, &tracekm1)); 81 PetscCall(ExpectedNumDofs_Total(1, order, 1, trimmed, 0, 1, &fiber1)); 82 } 83 *nDofs = (tracek * fiber0 + tracekm1 * fiber1) * nCopies; 84 } 85 PetscFunctionReturn(PETSC_SUCCESS); 86 } 87 88 static PetscErrorCode ExpectedNumDofs_Interior(PetscInt dim, PetscInt order, PetscInt formDegree, PetscBool trimmed, PetscInt tensor, PetscInt nCopies, PetscInt *nDofs) 89 { 90 PetscFunctionBegin; 91 formDegree = PetscAbsInt(formDegree); 92 /* see femtable.org for the source of most of these values */ 93 *nDofs = -1; 94 if (tensor == 0) { /* simplex spaces */ 95 if (trimmed) { 96 if (order + formDegree > dim) { 97 PetscInt eorder = order + formDegree - dim - 1; 98 PetscInt eformDegree = dim - formDegree; 99 PetscInt rnchooserk; 100 PetscInt rkchoosek; 101 102 PetscCall(PetscDTBinomialInt(eorder + dim, eorder + eformDegree, &rnchooserk)); 103 PetscCall(PetscDTBinomialInt(eorder + eformDegree, eformDegree, &rkchoosek)); 104 *nDofs = rnchooserk * rkchoosek * nCopies; 105 } else { 106 *nDofs = 0; 107 } 108 109 } else { 110 if (order + formDegree > dim) { 111 PetscInt eorder = order + formDegree - dim; 112 PetscInt eformDegree = dim - formDegree; 113 PetscInt rnchooserk; 114 PetscInt rkm1choosek; 115 116 PetscCall(PetscDTBinomialInt(eorder + dim, eorder + eformDegree, &rnchooserk)); 117 PetscCall(PetscDTBinomialInt(eorder + eformDegree - 1, eformDegree, &rkm1choosek)); 118 *nDofs = rnchooserk * rkm1choosek * nCopies; 119 } else { 120 *nDofs = 0; 121 } 122 } 123 } else if (tensor == 1) { /* hypercubes */ 124 if (dim < 2) { 125 PetscCall(ExpectedNumDofs_Interior(dim, order, formDegree, trimmed, 0, nCopies, nDofs)); 126 } else { 127 PetscInt tracek = 0; 128 PetscInt tracekm1 = 0; 129 PetscInt fiber0 = 0; 130 PetscInt fiber1 = 0; 131 132 if (formDegree < dim) { 133 PetscCall(ExpectedNumDofs_Interior(dim - 1, order, formDegree, trimmed, dim > 2 ? 1 : 0, 1, &tracek)); 134 PetscCall(ExpectedNumDofs_Interior(1, order, 0, trimmed, 0, 1, &fiber0)); 135 } 136 if (formDegree > 0) { 137 PetscCall(ExpectedNumDofs_Interior(dim - 1, order, formDegree - 1, trimmed, dim > 2 ? 1 : 0, 1, &tracekm1)); 138 PetscCall(ExpectedNumDofs_Interior(1, order, 1, trimmed, 0, 1, &fiber1)); 139 } 140 *nDofs = (tracek * fiber0 + tracekm1 * fiber1) * nCopies; 141 } 142 } else { /* prism */ 143 PetscInt tracek = 0; 144 PetscInt tracekm1 = 0; 145 PetscInt fiber0 = 0; 146 PetscInt fiber1 = 0; 147 148 if (formDegree < dim) { 149 PetscCall(ExpectedNumDofs_Interior(dim - 1, order, formDegree, trimmed, 0, 1, &tracek)); 150 PetscCall(ExpectedNumDofs_Interior(1, order, 0, trimmed, 0, 1, &fiber0)); 151 } 152 if (formDegree > 0) { 153 PetscCall(ExpectedNumDofs_Interior(dim - 1, order, formDegree - 1, trimmed, 0, 1, &tracekm1)); 154 PetscCall(ExpectedNumDofs_Interior(1, order, 1, trimmed, 0, 1, &fiber1)); 155 } 156 *nDofs = (tracek * fiber0 + tracekm1 * fiber1) * nCopies; 157 } 158 PetscFunctionReturn(PETSC_SUCCESS); 159 } 160 161 PetscErrorCode testLagrange(PetscHashLag lagTable, DM K, PetscInt dim, PetscInt order, PetscInt formDegree, PetscBool trimmed, PetscInt tensor, PetscBool continuous, PetscInt nCopies) 162 { 163 PetscDualSpace sp; 164 MPI_Comm comm = PETSC_COMM_SELF; 165 PetscInt Nk; 166 PetscHashLagKey key; 167 PetscHashIter iter; 168 PetscBool missing; 169 PetscInt spdim, spintdim, exspdim, exspintdim; 170 171 PetscFunctionBegin; 172 PetscCall(PetscDTBinomialInt(dim, PetscAbsInt(formDegree), &Nk)); 173 PetscCall(PetscDualSpaceCreate(comm, &sp)); 174 PetscCall(PetscDualSpaceSetType(sp, PETSCDUALSPACELAGRANGE)); 175 PetscCall(PetscDualSpaceSetDM(sp, K)); 176 PetscCall(PetscDualSpaceSetOrder(sp, order)); 177 PetscCall(PetscDualSpaceSetFormDegree(sp, formDegree)); 178 PetscCall(PetscDualSpaceSetNumComponents(sp, nCopies * Nk)); 179 PetscCall(PetscDualSpaceLagrangeSetContinuity(sp, continuous)); 180 PetscCall(PetscDualSpaceLagrangeSetTensor(sp, (PetscBool)tensor)); 181 PetscCall(PetscDualSpaceLagrangeSetTrimmed(sp, trimmed)); 182 PetscCall(PetscInfo(NULL, "Input: dim %" PetscInt_FMT ", order %" PetscInt_FMT ", trimmed %" PetscInt_FMT ", tensor %" PetscInt_FMT ", continuous %" PetscInt_FMT ", formDegree %" PetscInt_FMT ", nCopies %" PetscInt_FMT "\n", dim, order, (PetscInt)trimmed, tensor, (PetscInt)continuous, formDegree, nCopies)); 183 PetscCall(ExpectedNumDofs_Total(dim, order, formDegree, trimmed, tensor, nCopies, &exspdim)); 184 if (continuous && dim > 0 && order > 0) { 185 PetscCall(ExpectedNumDofs_Interior(dim, order, formDegree, trimmed, tensor, nCopies, &exspintdim)); 186 } else { 187 exspintdim = exspdim; 188 } 189 PetscCall(PetscDualSpaceSetUp(sp)); 190 PetscCall(PetscDualSpaceGetDimension(sp, &spdim)); 191 PetscCall(PetscDualSpaceGetInteriorDimension(sp, &spintdim)); 192 PetscCheck(spdim == exspdim, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Expected dual space dimension %" PetscInt_FMT ", got %" PetscInt_FMT, exspdim, spdim); 193 PetscCheck(spintdim == exspintdim, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Expected dual space interior dimension %" PetscInt_FMT ", got %" PetscInt_FMT, exspintdim, spintdim); 194 key.dim = dim; 195 key.formDegree = formDegree; 196 PetscCall(PetscDualSpaceGetOrder(sp, &key.order)); 197 PetscCall(PetscDualSpaceLagrangeGetContinuity(sp, &key.continuous)); 198 if (tensor == 2) { 199 key.tensor = 2; 200 } else { 201 PetscBool bTensor; 202 203 PetscCall(PetscDualSpaceLagrangeGetTensor(sp, &bTensor)); 204 key.tensor = bTensor; 205 } 206 PetscCall(PetscDualSpaceLagrangeGetTrimmed(sp, &key.trimmed)); 207 PetscCall(PetscInfo(NULL, "After setup: order %" PetscInt_FMT ", trimmed %" PetscInt_FMT ", tensor %" PetscInt_FMT ", continuous %" PetscInt_FMT "\n", key.order, (PetscInt)key.trimmed, key.tensor, (PetscInt)key.continuous)); 208 PetscCall(PetscHashLagPut(lagTable, key, &iter, &missing)); 209 if (missing) { 210 DMPolytopeType type; 211 212 PetscCall(DMPlexGetCellType(K, 0, &type)); 213 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_SELF, "New space: %s, order %" PetscInt_FMT ", trimmed %" PetscInt_FMT ", tensor %" PetscInt_FMT ", continuous %" PetscInt_FMT ", form degree %" PetscInt_FMT "\n", DMPolytopeTypes[type], order, (PetscInt)trimmed, tensor, (PetscInt)continuous, formDegree)); 214 PetscCall(PetscViewerASCIIPushTab(PETSC_VIEWER_STDOUT_SELF)); 215 { 216 PetscQuadrature intNodes, allNodes; 217 Mat intMat, allMat; 218 MatInfo info; 219 PetscInt i, j, nodeIdxDim, nodeVecDim, nNodes; 220 const PetscInt *nodeIdx; 221 const PetscReal *nodeVec; 222 223 PetscDualSpace_Lag *lag = (PetscDualSpace_Lag *)sp->data; 224 225 PetscCall(PetscLagNodeIndicesGetData_Internal(lag->allNodeIndices, &nodeIdxDim, &nodeVecDim, &nNodes, &nodeIdx, &nodeVec)); 226 PetscCheck(nodeVecDim == Nk, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Incorrect nodeVecDim"); 227 PetscCheck(nNodes == spdim, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Incorrect nNodes"); 228 229 PetscCall(PetscDualSpaceGetAllData(sp, &allNodes, &allMat)); 230 231 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_SELF, "All nodes:\n")); 232 PetscCall(PetscViewerASCIIPushTab(PETSC_VIEWER_STDOUT_SELF)); 233 PetscCall(PetscQuadratureView(allNodes, PETSC_VIEWER_STDOUT_SELF)); 234 PetscCall(PetscViewerASCIIPopTab(PETSC_VIEWER_STDOUT_SELF)); 235 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_SELF, "All node indices:\n")); 236 for (i = 0; i < spdim; i++) { 237 PetscCall(PetscPrintf(PETSC_COMM_SELF, "(")); 238 for (j = 0; j < nodeIdxDim; j++) PetscCall(PetscPrintf(PETSC_COMM_SELF, " %" PetscInt_FMT ",", nodeIdx[i * nodeIdxDim + j])); 239 PetscCall(PetscPrintf(PETSC_COMM_SELF, "): [")); 240 for (j = 0; j < nodeVecDim; j++) PetscCall(PetscPrintf(PETSC_COMM_SELF, " %g,", (double)nodeVec[i * nodeVecDim + j])); 241 PetscCall(PetscPrintf(PETSC_COMM_SELF, "]\n")); 242 } 243 244 PetscCall(MatGetInfo(allMat, MAT_LOCAL, &info)); 245 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_SELF, "All matrix: %" PetscInt_FMT " nonzeros\n", (PetscInt)info.nz_used)); 246 247 PetscCall(PetscDualSpaceGetInteriorData(sp, &intNodes, &intMat)); 248 if (intMat && intMat != allMat) { 249 PetscInt intNodeIdxDim, intNodeVecDim, intNnodes; 250 const PetscInt *intNodeIdx; 251 const PetscReal *intNodeVec; 252 PetscBool same; 253 254 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_SELF, "Interior nodes:\n")); 255 PetscCall(PetscViewerASCIIPushTab(PETSC_VIEWER_STDOUT_SELF)); 256 PetscCall(PetscQuadratureView(intNodes, PETSC_VIEWER_STDOUT_SELF)); 257 PetscCall(PetscViewerASCIIPopTab(PETSC_VIEWER_STDOUT_SELF)); 258 259 PetscCall(MatGetInfo(intMat, MAT_LOCAL, &info)); 260 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_SELF, "Interior matrix: %" PetscInt_FMT " nonzeros\n", (PetscInt)info.nz_used)); 261 PetscCall(PetscLagNodeIndicesGetData_Internal(lag->intNodeIndices, &intNodeIdxDim, &intNodeVecDim, &intNnodes, &intNodeIdx, &intNodeVec)); 262 PetscCheck(intNodeIdxDim == nodeIdxDim && intNodeVecDim == nodeVecDim, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Interior node indices not the same shale as all node indices"); 263 PetscCheck(intNnodes == spintdim, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Incorrect interior nNodes"); 264 PetscCall(PetscArraycmp(intNodeIdx, nodeIdx, nodeIdxDim * intNnodes, &same)); 265 PetscCheck(same, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Interior node indices not the same as start of all node indices"); 266 PetscCall(PetscArraycmp(intNodeVec, nodeVec, nodeVecDim * intNnodes, &same)); 267 PetscCheck(same, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Interior node vectors not the same as start of all node vectors"); 268 } else if (intMat) { 269 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_SELF, "Interior data is the same as all data\n")); 270 PetscCheck(intNodes == allNodes, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Interior nodes should be the same as all nodes"); 271 PetscCheck(lag->intNodeIndices == lag->allNodeIndices, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Interior node indices should be the same as all node indices"); 272 } 273 } 274 if (dim <= 2 && spintdim) { 275 PetscInt numFaces, o; 276 277 { 278 DMPolytopeType ct; 279 /* The number of arrangements is no longer based on the number of faces */ 280 PetscCall(DMPlexGetCellType(K, 0, &ct)); 281 numFaces = DMPolytopeTypeGetNumArrangments(ct) / 2; 282 } 283 for (o = -numFaces; o < numFaces; ++o) { 284 Mat symMat; 285 286 PetscCall(PetscDualSpaceCreateInteriorSymmetryMatrix_Lagrange(sp, o, &symMat)); 287 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_SELF, "Interior node symmetry matrix for orientation %" PetscInt_FMT ":\n", o)); 288 PetscCall(PetscViewerASCIIPushTab(PETSC_VIEWER_STDOUT_SELF)); 289 PetscCall(MatView(symMat, PETSC_VIEWER_STDOUT_SELF)); 290 PetscCall(PetscViewerASCIIPopTab(PETSC_VIEWER_STDOUT_SELF)); 291 PetscCall(MatDestroy(&symMat)); 292 } 293 } 294 PetscCall(PetscViewerASCIIPopTab(PETSC_VIEWER_STDOUT_SELF)); 295 } 296 PetscCall(PetscDualSpaceDestroy(&sp)); 297 PetscFunctionReturn(PETSC_SUCCESS); 298 } 299 300 int main(int argc, char **argv) 301 { 302 PetscInt dim; 303 PetscHashLag lagTable; 304 PetscInt tensorCell; 305 PetscInt order, ordermin, ordermax; 306 PetscBool continuous; 307 PetscBool trimmed; 308 DM dm; 309 310 PetscFunctionBeginUser; 311 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 312 dim = 3; 313 tensorCell = 0; 314 continuous = PETSC_FALSE; 315 trimmed = PETSC_FALSE; 316 PetscOptionsBegin(PETSC_COMM_WORLD, "", "Options for PETSCDUALSPACELAGRANGE test", "none"); 317 PetscCall(PetscOptionsRangeInt("-dim", "The spatial dimension", "ex1.c", dim, &dim, NULL, 0, 3)); 318 PetscCall(PetscOptionsRangeInt("-tensor", "(0) simplex (1) hypercube (2) wedge", "ex1.c", tensorCell, &tensorCell, NULL, 0, 2)); 319 PetscCall(PetscOptionsBool("-continuous", "Whether the dual space has continuity", "ex1.c", continuous, &continuous, NULL)); 320 PetscCall(PetscOptionsBool("-trimmed", "Whether the dual space matches a trimmed polynomial space", "ex1.c", trimmed, &trimmed, NULL)); 321 PetscOptionsEnd(); 322 PetscCall(PetscHashLagCreate(&lagTable)); 323 324 if (tensorCell < 2) { 325 PetscCall(DMPlexCreateReferenceCell(PETSC_COMM_SELF, DMPolytopeTypeSimpleShape(dim, (PetscBool)!tensorCell), &dm)); 326 } else { 327 PetscCall(DMPlexCreateReferenceCell(PETSC_COMM_SELF, DM_POLYTOPE_TRI_PRISM, &dm)); 328 } 329 ordermin = trimmed ? 1 : 0; 330 ordermax = tensorCell == 2 ? 4 : tensorCell == 1 ? 3 : dim + 2; 331 for (order = ordermin; order <= ordermax; order++) { 332 PetscInt formDegree; 333 334 for (formDegree = PetscMin(0, -dim + 1); formDegree <= dim; formDegree++) { 335 PetscInt nCopies; 336 337 for (nCopies = 1; nCopies <= 3; nCopies++) PetscCall(testLagrange(lagTable, dm, dim, order, formDegree, trimmed, (PetscBool)tensorCell, continuous, nCopies)); 338 } 339 } 340 PetscCall(DMDestroy(&dm)); 341 PetscCall(PetscHashLagDestroy(&lagTable)); 342 PetscCall(PetscFinalize()); 343 return 0; 344 } 345 346 /*TEST 347 348 test: 349 suffix: 0 350 args: -dim 0 351 352 test: 353 suffix: 1_discontinuous_full 354 args: -dim 1 -continuous 0 -trimmed 0 355 356 test: 357 suffix: 1_continuous_full 358 args: -dim 1 -continuous 1 -trimmed 0 359 360 test: 361 suffix: 2_simplex_discontinuous_full 362 args: -dim 2 -tensor 0 -continuous 0 -trimmed 0 363 364 test: 365 suffix: 2_simplex_continuous_full 366 args: -dim 2 -tensor 0 -continuous 1 -trimmed 0 367 368 test: 369 suffix: 2_tensor_discontinuous_full 370 args: -dim 2 -tensor 1 -continuous 0 -trimmed 0 371 372 test: 373 suffix: 2_tensor_continuous_full 374 args: -dim 2 -tensor 1 -continuous 1 -trimmed 0 375 376 test: 377 suffix: 3_simplex_discontinuous_full 378 args: -dim 3 -tensor 0 -continuous 0 -trimmed 0 379 380 test: 381 suffix: 3_simplex_continuous_full 382 args: -dim 3 -tensor 0 -continuous 1 -trimmed 0 383 384 test: 385 suffix: 3_tensor_discontinuous_full 386 args: -dim 3 -tensor 1 -continuous 0 -trimmed 0 387 388 test: 389 suffix: 3_tensor_continuous_full 390 args: -dim 3 -tensor 1 -continuous 1 -trimmed 0 391 392 test: 393 suffix: 3_wedge_discontinuous_full 394 args: -dim 3 -tensor 2 -continuous 0 -trimmed 0 395 396 test: 397 suffix: 3_wedge_continuous_full 398 args: -dim 3 -tensor 2 -continuous 1 -trimmed 0 399 400 test: 401 suffix: 1_discontinuous_trimmed 402 args: -dim 1 -continuous 0 -trimmed 1 403 404 test: 405 suffix: 1_continuous_trimmed 406 args: -dim 1 -continuous 1 -trimmed 1 407 408 test: 409 suffix: 2_simplex_discontinuous_trimmed 410 args: -dim 2 -tensor 0 -continuous 0 -trimmed 1 411 412 test: 413 suffix: 2_simplex_continuous_trimmed 414 args: -dim 2 -tensor 0 -continuous 1 -trimmed 1 415 416 test: 417 suffix: 2_tensor_discontinuous_trimmed 418 args: -dim 2 -tensor 1 -continuous 0 -trimmed 1 419 420 test: 421 suffix: 2_tensor_continuous_trimmed 422 args: -dim 2 -tensor 1 -continuous 1 -trimmed 1 423 424 test: 425 suffix: 3_simplex_discontinuous_trimmed 426 args: -dim 3 -tensor 0 -continuous 0 -trimmed 1 427 428 test: 429 suffix: 3_simplex_continuous_trimmed 430 args: -dim 3 -tensor 0 -continuous 1 -trimmed 1 431 432 test: 433 suffix: 3_tensor_discontinuous_trimmed 434 args: -dim 3 -tensor 1 -continuous 0 -trimmed 1 435 436 test: 437 suffix: 3_tensor_continuous_trimmed 438 args: -dim 3 -tensor 1 -continuous 1 -trimmed 1 439 440 test: 441 suffix: 3_wedge_discontinuous_trimmed 442 args: -dim 3 -tensor 2 -continuous 0 -trimmed 1 443 444 test: 445 suffix: 3_wedge_continuous_trimmed 446 args: -dim 3 -tensor 2 -continuous 1 -trimmed 1 447 448 TEST*/ 449