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