1 static char help[] = "Tests for DMLabel\n\n"; 2 3 #include <petscdmplex.h> 4 #include <petsc/private/dmimpl.h> 5 6 static PetscErrorCode TestInsertion() 7 { 8 DMLabel label, label2; 9 const PetscInt values[5] = {0, 3, 4, -1, 176}, N = 10000; 10 PetscInt i, v; 11 12 PetscFunctionBegin; 13 PetscCall(DMLabelCreate(PETSC_COMM_SELF, "Test Label", &label)); 14 PetscCall(DMLabelSetDefaultValue(label, -100)); 15 for (i = 0; i < N; ++i) { 16 PetscCall(DMLabelSetValue(label, i, values[i%5])); 17 } 18 /* Test get in hash mode */ 19 for (i = 0; i < N; ++i) { 20 PetscInt val; 21 22 PetscCall(DMLabelGetValue(label, i, &val)); 23 PetscCheck(val == values[i%5],PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Value %" PetscInt_FMT " for point %" PetscInt_FMT " should be %" PetscInt_FMT, val, i, values[i%5]); 24 } 25 /* Test stratum */ 26 for (v = 0; v < 5; ++v) { 27 IS stratum; 28 const PetscInt *points; 29 PetscInt n; 30 31 PetscCall(DMLabelGetStratumIS(label, values[v], &stratum)); 32 PetscCheck(stratum,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Stratum %" PetscInt_FMT " is empty!", v); 33 PetscCall(ISGetIndices(stratum, &points)); 34 PetscCall(ISGetLocalSize(stratum, &n)); 35 for (i = 0; i < n; ++i) { 36 PetscCheck(points[i] == i*5+v,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Point %" PetscInt_FMT " should be %" PetscInt_FMT, points[i], i*5+v); 37 } 38 PetscCall(ISRestoreIndices(stratum, &points)); 39 PetscCall(ISDestroy(&stratum)); 40 } 41 /* Test get in array mode */ 42 for (i = 0; i < N; ++i) { 43 PetscInt val; 44 45 PetscCall(DMLabelGetValue(label, i, &val)); 46 PetscCheck(val == values[i%5],PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Value %" PetscInt_FMT " should be %" PetscInt_FMT, val, values[i%5]); 47 } 48 /* Test Duplicate */ 49 PetscCall(DMLabelDuplicate(label, &label2)); 50 for (i = 0; i < N; ++i) { 51 PetscInt val; 52 53 PetscCall(DMLabelGetValue(label2, i, &val)); 54 PetscCheck(val == values[i%5],PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Value %" PetscInt_FMT " should be %" PetscInt_FMT, val, values[i%5]); 55 } 56 PetscCall(DMLabelDestroy(&label2)); 57 PetscCall(DMLabelDestroy(&label)); 58 PetscFunctionReturn(0); 59 } 60 61 static PetscErrorCode TestEmptyStrata(MPI_Comm comm) 62 { 63 DM dm, dmDist; 64 PetscPartitioner part; 65 PetscInt c0[6] = {2,3,6,7,9,11}; 66 PetscInt c1[6] = {4,5,7,8,10,12}; 67 PetscInt c2[4] = {13,15,19,21}; 68 PetscInt c3[4] = {14,16,20,22}; 69 PetscInt c4[4] = {15,17,21,23}; 70 PetscInt c5[4] = {16,18,22,24}; 71 PetscInt c6[4] = {13,14,19,20}; 72 PetscInt c7[4] = {15,16,21,22}; 73 PetscInt c8[4] = {17,18,23,24}; 74 PetscInt c9[4] = {13,14,15,16}; 75 PetscInt c10[4] = {15,16,17,18}; 76 PetscInt c11[4] = {19,20,21,22}; 77 PetscInt c12[4] = {21,22,23,24}; 78 PetscInt dim = 3; 79 PetscMPIInt rank; 80 81 PetscFunctionBegin; 82 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 83 /* A 3D box with two adjacent cells, sharing one face and four vertices */ 84 PetscCall(DMCreate(comm, &dm)); 85 PetscCall(DMSetType(dm, DMPLEX)); 86 PetscCall(DMSetDimension(dm, dim)); 87 if (rank == 0) { 88 PetscCall(DMPlexSetChart(dm, 0, 25)); 89 PetscCall(DMPlexSetConeSize(dm, 0, 6)); 90 PetscCall(DMPlexSetConeSize(dm, 1, 6)); 91 PetscCall(DMPlexSetConeSize(dm, 2, 4)); 92 PetscCall(DMPlexSetConeSize(dm, 3, 4)); 93 PetscCall(DMPlexSetConeSize(dm, 4, 4)); 94 PetscCall(DMPlexSetConeSize(dm, 5, 4)); 95 PetscCall(DMPlexSetConeSize(dm, 6, 4)); 96 PetscCall(DMPlexSetConeSize(dm, 7, 4)); 97 PetscCall(DMPlexSetConeSize(dm, 8, 4)); 98 PetscCall(DMPlexSetConeSize(dm, 9, 4)); 99 PetscCall(DMPlexSetConeSize(dm, 10, 4)); 100 PetscCall(DMPlexSetConeSize(dm, 11, 4)); 101 PetscCall(DMPlexSetConeSize(dm, 12, 4)); 102 } 103 PetscCall(DMSetUp(dm)); 104 if (rank == 0) { 105 PetscCall(DMPlexSetCone(dm, 0, c0)); 106 PetscCall(DMPlexSetCone(dm, 1, c1)); 107 PetscCall(DMPlexSetCone(dm, 2, c2)); 108 PetscCall(DMPlexSetCone(dm, 3, c3)); 109 PetscCall(DMPlexSetCone(dm, 4, c4)); 110 PetscCall(DMPlexSetCone(dm, 5, c5)); 111 PetscCall(DMPlexSetCone(dm, 6, c6)); 112 PetscCall(DMPlexSetCone(dm, 7, c7)); 113 PetscCall(DMPlexSetCone(dm, 8, c8)); 114 PetscCall(DMPlexSetCone(dm, 9, c9)); 115 PetscCall(DMPlexSetCone(dm, 10, c10)); 116 PetscCall(DMPlexSetCone(dm, 11, c11)); 117 PetscCall(DMPlexSetCone(dm, 12, c12)); 118 } 119 PetscCall(DMPlexSymmetrize(dm)); 120 /* Create a user managed depth label, so that we can leave out edges */ 121 { 122 DMLabel label; 123 PetscInt numValues, maxValues = 0, v; 124 125 PetscCall(DMCreateLabel(dm, "depth")); 126 PetscCall(DMPlexGetDepthLabel(dm, &label)); 127 if (rank == 0) { 128 PetscInt i; 129 130 for (i = 0; i < 25; ++i) { 131 if (i < 2) PetscCall(DMLabelSetValue(label, i, 3)); 132 else if (i < 13) PetscCall(DMLabelSetValue(label, i, 2)); 133 else { 134 if (i==13) PetscCall(DMLabelAddStratum(label, 1)); 135 PetscCall(DMLabelSetValue(label, i, 0)); 136 } 137 } 138 } 139 PetscCall(DMLabelGetNumValues(label, &numValues)); 140 PetscCallMPI(MPI_Allreduce(&numValues, &maxValues, 1, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject) dm))); 141 for (v = numValues; v < maxValues; ++v) PetscCall(DMLabelAddStratum(label,v)); 142 } 143 { 144 DMLabel label; 145 PetscCall(DMPlexGetDepthLabel(dm, &label)); 146 PetscCall(DMLabelView(label, PETSC_VIEWER_STDOUT_(comm))); 147 } 148 PetscCall(DMPlexGetPartitioner(dm,&part)); 149 PetscCall(PetscPartitionerSetFromOptions(part)); 150 PetscCall(DMPlexDistribute(dm, 1, NULL, &dmDist)); 151 if (dmDist) { 152 PetscCall(DMDestroy(&dm)); 153 dm = dmDist; 154 } 155 { 156 DMLabel label; 157 PetscCall(DMPlexGetDepthLabel(dm, &label)); 158 PetscCall(DMLabelView(label, PETSC_VIEWER_STDOUT_(comm))); 159 } 160 /* Create a cell vector */ 161 { 162 Vec v; 163 PetscSection s; 164 PetscInt numComp[] = {1}; 165 PetscInt dof[] = {0,0,0,1}; 166 PetscInt N; 167 168 PetscCall(DMSetNumFields(dm, 1)); 169 PetscCall(DMPlexCreateSection(dm, NULL, numComp, dof, 0, NULL, NULL, NULL, NULL, &s)); 170 PetscCall(DMSetLocalSection(dm, s)); 171 PetscCall(PetscSectionDestroy(&s)); 172 PetscCall(DMCreateGlobalVector(dm, &v)); 173 PetscCall(VecGetSize(v, &N)); 174 if (N != 2) { 175 PetscCall(DMView(dm, PETSC_VIEWER_STDOUT_(comm))); 176 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "FAIL: Vector size %" PetscInt_FMT " != 2", N); 177 } 178 PetscCall(VecDestroy(&v)); 179 } 180 PetscCall(DMDestroy(&dm)); 181 PetscFunctionReturn(0); 182 } 183 184 static PetscErrorCode TestDistribution(MPI_Comm comm) 185 { 186 DM dm, dmDist; 187 PetscPartitioner part; 188 DMLabel label; 189 char filename[PETSC_MAX_PATH_LEN]; 190 const char *name = "test label"; 191 PetscInt overlap = 0, cStart, cEnd, c; 192 PetscMPIInt rank; 193 PetscBool flg; 194 195 PetscFunctionBegin; 196 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 197 PetscCall(PetscOptionsGetString(NULL, NULL, "-filename", filename, sizeof(filename), &flg)); 198 if (!flg) PetscFunctionReturn(0); 199 PetscCall(PetscOptionsGetInt(NULL, NULL, "-overlap", &overlap, NULL)); 200 PetscCall(DMPlexCreateFromFile(comm, filename, "ex11_plex", PETSC_TRUE, &dm)); 201 PetscCall(DMSetBasicAdjacency(dm, PETSC_TRUE, PETSC_FALSE)); 202 PetscCall(DMCreateLabel(dm, name)); 203 PetscCall(DMGetLabel(dm, name, &label)); 204 PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd)); 205 for (c = cStart; c < cEnd; ++c) { 206 PetscCall(DMLabelSetValue(label, c, c)); 207 } 208 PetscCall(DMLabelView(label, PETSC_VIEWER_STDOUT_WORLD)); 209 PetscCall(DMPlexGetPartitioner(dm,&part)); 210 PetscCall(PetscPartitionerSetFromOptions(part)); 211 PetscCall(DMPlexDistribute(dm, overlap, NULL, &dmDist)); 212 if (dmDist) { 213 PetscCall(DMDestroy(&dm)); 214 dm = dmDist; 215 } 216 PetscCall(PetscObjectSetName((PetscObject) dm, "Mesh")); 217 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view")); 218 PetscCall(DMGetLabel(dm, name, &label)); 219 PetscCall(DMLabelView(label, PETSC_VIEWER_STDOUT_WORLD)); 220 PetscCall(DMDestroy(&dm)); 221 PetscFunctionReturn(0); 222 } 223 224 static PetscErrorCode TestUniversalLabel(MPI_Comm comm) 225 { 226 DM dm1, dm2; 227 DMLabel bd1, bd2, ulabel; 228 DMUniversalLabel universal; 229 PetscInt pStart, pEnd, p; 230 PetscBool run = PETSC_FALSE, notFile; 231 232 PetscFunctionBeginUser; 233 PetscCall(PetscOptionsGetBool(NULL, NULL, "-universal", &run, NULL)); 234 if (!run) PetscFunctionReturn(0); 235 236 char filename[PETSC_MAX_PATH_LEN]; 237 PetscBool flg; 238 239 PetscCall(PetscOptionsGetString(NULL, NULL, "-filename", filename, sizeof(filename), &flg)); 240 if (flg) { 241 PetscCall(DMPlexCreateFromFile(comm, filename, "ex11_plex", PETSC_TRUE, &dm1)); 242 } else { 243 PetscCall(DMCreate(comm, &dm1)); 244 PetscCall(DMSetType(dm1, DMPLEX)); 245 PetscCall(DMSetFromOptions(dm1)); 246 } 247 PetscCall(DMHasLabel(dm1, "marker", ¬File)); 248 if (notFile) { 249 PetscCall(DMCreateLabel(dm1, "Boundary Faces")); 250 PetscCall(DMGetLabel(dm1, "Boundary Faces", &bd1)); 251 PetscCall(DMPlexMarkBoundaryFaces(dm1, 13, bd1)); 252 PetscCall(DMCreateLabel(dm1, "Boundary")); 253 PetscCall(DMGetLabel(dm1, "Boundary", &bd2)); 254 PetscCall(DMPlexMarkBoundaryFaces(dm1, 121, bd2)); 255 PetscCall(DMPlexLabelComplete(dm1, bd2)); 256 } 257 PetscCall(PetscObjectSetName((PetscObject) dm1, "First Mesh")); 258 PetscCall(DMViewFromOptions(dm1, NULL, "-dm_view")); 259 260 PetscCall(DMUniversalLabelCreate(dm1, &universal)); 261 PetscCall(DMUniversalLabelGetLabel(universal, &ulabel)); 262 PetscCall(PetscObjectViewFromOptions((PetscObject) ulabel, NULL, "-universal_view")); 263 264 if (!notFile) { 265 PetscInt Nl, l; 266 267 PetscCall(DMClone(dm1, &dm2)); 268 PetscCall(DMGetNumLabels(dm2, &Nl)); 269 for (l = Nl-1; l >= 0; --l) { 270 PetscBool isdepth, iscelltype; 271 const char *name; 272 273 PetscCall(DMGetLabelName(dm2, l, &name)); 274 PetscCall(PetscStrncmp(name, "depth", 6, &isdepth)); 275 PetscCall(PetscStrncmp(name, "celltype", 9, &iscelltype)); 276 if (!isdepth && !iscelltype) PetscCall(DMRemoveLabel(dm2, name, NULL)); 277 } 278 } else { 279 PetscCall(DMCreate(comm, &dm2)); 280 PetscCall(DMSetType(dm2, DMPLEX)); 281 PetscCall(DMSetFromOptions(dm2)); 282 } 283 PetscCall(PetscObjectSetName((PetscObject) dm2, "Second Mesh")); 284 PetscCall(DMUniversalLabelCreateLabels(universal, PETSC_TRUE, dm2)); 285 PetscCall(DMPlexGetChart(dm2, &pStart, &pEnd)); 286 for (p = pStart; p < pEnd; ++p) { 287 PetscInt val; 288 289 PetscCall(DMLabelGetValue(ulabel, p, &val)); 290 if (val < 0) continue; 291 PetscCall(DMUniversalLabelSetLabelValue(universal, dm2, PETSC_TRUE, p, val)); 292 } 293 PetscCall(DMViewFromOptions(dm2, NULL, "-dm_view")); 294 295 PetscCall(DMUniversalLabelDestroy(&universal)); 296 PetscCall(DMDestroy(&dm1)); 297 PetscCall(DMDestroy(&dm2)); 298 PetscFunctionReturn(0); 299 } 300 301 int main(int argc, char **argv) 302 { 303 304 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 305 /*PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));*/ 306 PetscCall(TestInsertion()); 307 PetscCall(TestEmptyStrata(PETSC_COMM_WORLD)); 308 PetscCall(TestDistribution(PETSC_COMM_WORLD)); 309 PetscCall(TestUniversalLabel(PETSC_COMM_WORLD)); 310 PetscCall(PetscFinalize()); 311 return 0; 312 } 313 314 /*TEST 315 316 test: 317 suffix: 0 318 requires: triangle 319 test: 320 suffix: 1 321 requires: triangle 322 nsize: 2 323 args: -petscpartitioner_type simple 324 325 testset: 326 suffix: gmsh 327 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.msh -petscpartitioner_type simple 328 test: 329 suffix: 1 330 nsize: 1 331 test: 332 suffix: 2 333 nsize: 2 334 335 testset: 336 suffix: exodusii 337 requires: exodusii 338 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/2Dgrd.exo -petscpartitioner_type simple 339 test: 340 suffix: 1 341 nsize: 1 342 test: 343 suffix: 2 344 nsize: 2 345 346 test: 347 suffix: univ 348 requires: triangle 349 args: -universal -dm_view -universal_view 350 351 test: 352 # Note that the labels differ because we have multiply-marked some points during EGADS creation 353 suffix: univ_egads_sphere 354 requires: egads 355 args: -universal -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/unit_sphere.egadslite -dm_view -universal_view 356 357 test: 358 # Note that the labels differ because we have multiply-marked some points during EGADS creation 359 suffix: univ_egads_ball 360 requires: egads ctetgen 361 args: -universal -dm_plex_boundary_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/unit_sphere.egadslite -dm_view -universal_view 362 363 TEST*/ 364