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