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 CHKERRQ(DMLabelCreate(PETSC_COMM_SELF, "Test Label", &label)); 14 CHKERRQ(DMLabelSetDefaultValue(label, -100)); 15 for (i = 0; i < N; ++i) { 16 CHKERRQ(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 CHKERRQ(DMLabelGetValue(label, i, &val)); 23 PetscCheckFalse(val != values[i%5],PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Value %d for point %d should be %d", 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 CHKERRQ(DMLabelGetStratumIS(label, values[v], &stratum)); 32 PetscCheck(stratum,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Stratum %d is empty!", v); 33 CHKERRQ(ISGetIndices(stratum, &points)); 34 CHKERRQ(ISGetLocalSize(stratum, &n)); 35 for (i = 0; i < n; ++i) { 36 PetscCheckFalse(points[i] != i*5+v,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Point %d should be %d", points[i], i*5+v); 37 } 38 CHKERRQ(ISRestoreIndices(stratum, &points)); 39 CHKERRQ(ISDestroy(&stratum)); 40 } 41 /* Test get in array mode */ 42 for (i = 0; i < N; ++i) { 43 PetscInt val; 44 45 CHKERRQ(DMLabelGetValue(label, i, &val)); 46 PetscCheckFalse(val != values[i%5],PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Value %d should be %d", val, values[i%5]); 47 } 48 /* Test Duplicate */ 49 CHKERRQ(DMLabelDuplicate(label, &label2)); 50 for (i = 0; i < N; ++i) { 51 PetscInt val; 52 53 CHKERRQ(DMLabelGetValue(label2, i, &val)); 54 PetscCheckFalse(val != values[i%5],PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Value %d should be %d", val, values[i%5]); 55 } 56 CHKERRQ(DMLabelDestroy(&label2)); 57 CHKERRQ(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 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 83 /* A 3D box with two adjacent cells, sharing one face and four vertices */ 84 CHKERRQ(DMCreate(comm, &dm)); 85 CHKERRQ(DMSetType(dm, DMPLEX)); 86 CHKERRQ(DMSetDimension(dm, dim)); 87 if (rank == 0) { 88 CHKERRQ(DMPlexSetChart(dm, 0, 25)); 89 CHKERRQ(DMPlexSetConeSize(dm, 0, 6)); 90 CHKERRQ(DMPlexSetConeSize(dm, 1, 6)); 91 CHKERRQ(DMPlexSetConeSize(dm, 2, 4)); 92 CHKERRQ(DMPlexSetConeSize(dm, 3, 4)); 93 CHKERRQ(DMPlexSetConeSize(dm, 4, 4)); 94 CHKERRQ(DMPlexSetConeSize(dm, 5, 4)); 95 CHKERRQ(DMPlexSetConeSize(dm, 6, 4)); 96 CHKERRQ(DMPlexSetConeSize(dm, 7, 4)); 97 CHKERRQ(DMPlexSetConeSize(dm, 8, 4)); 98 CHKERRQ(DMPlexSetConeSize(dm, 9, 4)); 99 CHKERRQ(DMPlexSetConeSize(dm, 10, 4)); 100 CHKERRQ(DMPlexSetConeSize(dm, 11, 4)); 101 CHKERRQ(DMPlexSetConeSize(dm, 12, 4)); 102 } 103 CHKERRQ(DMSetUp(dm)); 104 if (rank == 0) { 105 CHKERRQ(DMPlexSetCone(dm, 0, c0)); 106 CHKERRQ(DMPlexSetCone(dm, 1, c1)); 107 CHKERRQ(DMPlexSetCone(dm, 2, c2)); 108 CHKERRQ(DMPlexSetCone(dm, 3, c3)); 109 CHKERRQ(DMPlexSetCone(dm, 4, c4)); 110 CHKERRQ(DMPlexSetCone(dm, 5, c5)); 111 CHKERRQ(DMPlexSetCone(dm, 6, c6)); 112 CHKERRQ(DMPlexSetCone(dm, 7, c7)); 113 CHKERRQ(DMPlexSetCone(dm, 8, c8)); 114 CHKERRQ(DMPlexSetCone(dm, 9, c9)); 115 CHKERRQ(DMPlexSetCone(dm, 10, c10)); 116 CHKERRQ(DMPlexSetCone(dm, 11, c11)); 117 CHKERRQ(DMPlexSetCone(dm, 12, c12)); 118 } 119 CHKERRQ(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 CHKERRQ(DMCreateLabel(dm, "depth")); 126 CHKERRQ(DMPlexGetDepthLabel(dm, &label)); 127 if (rank == 0) { 128 PetscInt i; 129 130 for (i = 0; i < 25; ++i) { 131 if (i < 2) CHKERRQ(DMLabelSetValue(label, i, 3)); 132 else if (i < 13) CHKERRQ(DMLabelSetValue(label, i, 2)); 133 else { 134 if (i==13) CHKERRQ(DMLabelAddStratum(label, 1)); 135 CHKERRQ(DMLabelSetValue(label, i, 0)); 136 } 137 } 138 } 139 CHKERRQ(DMLabelGetNumValues(label, &numValues)); 140 CHKERRMPI(MPI_Allreduce(&numValues, &maxValues, 1, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject) dm))); 141 for (v = numValues; v < maxValues; ++v) CHKERRQ(DMLabelAddStratum(label,v)); 142 } 143 { 144 DMLabel label; 145 CHKERRQ(DMPlexGetDepthLabel(dm, &label)); 146 CHKERRQ(DMLabelView(label, PETSC_VIEWER_STDOUT_(comm))); 147 } 148 CHKERRQ(DMPlexGetPartitioner(dm,&part)); 149 CHKERRQ(PetscPartitionerSetFromOptions(part)); 150 CHKERRQ(DMPlexDistribute(dm, 1, NULL, &dmDist)); 151 if (dmDist) { 152 CHKERRQ(DMDestroy(&dm)); 153 dm = dmDist; 154 } 155 { 156 DMLabel label; 157 CHKERRQ(DMPlexGetDepthLabel(dm, &label)); 158 CHKERRQ(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 CHKERRQ(DMSetNumFields(dm, 1)); 169 CHKERRQ(DMPlexCreateSection(dm, NULL, numComp, dof, 0, NULL, NULL, NULL, NULL, &s)); 170 CHKERRQ(DMSetLocalSection(dm, s)); 171 CHKERRQ(PetscSectionDestroy(&s)); 172 CHKERRQ(DMCreateGlobalVector(dm, &v)); 173 CHKERRQ(VecGetSize(v, &N)); 174 if (N != 2) { 175 CHKERRQ(DMView(dm, PETSC_VIEWER_STDOUT_(comm))); 176 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "FAIL: Vector size %d != 2", N); 177 } 178 CHKERRQ(VecDestroy(&v)); 179 } 180 CHKERRQ(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 CHKERRMPI(MPI_Comm_rank(comm, &rank)); 197 CHKERRQ(PetscOptionsGetString(NULL, NULL, "-filename", filename, sizeof(filename), &flg)); 198 if (!flg) PetscFunctionReturn(0); 199 CHKERRQ(PetscOptionsGetInt(NULL, NULL, "-overlap", &overlap, NULL)); 200 CHKERRQ(DMPlexCreateFromFile(comm, filename, "ex11_plex", PETSC_TRUE, &dm)); 201 CHKERRQ(DMSetBasicAdjacency(dm, PETSC_TRUE, PETSC_FALSE)); 202 CHKERRQ(DMCreateLabel(dm, name)); 203 CHKERRQ(DMGetLabel(dm, name, &label)); 204 CHKERRQ(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd)); 205 for (c = cStart; c < cEnd; ++c) { 206 CHKERRQ(DMLabelSetValue(label, c, c)); 207 } 208 CHKERRQ(DMLabelView(label, PETSC_VIEWER_STDOUT_WORLD)); 209 CHKERRQ(DMPlexGetPartitioner(dm,&part)); 210 CHKERRQ(PetscPartitionerSetFromOptions(part)); 211 CHKERRQ(DMPlexDistribute(dm, overlap, NULL, &dmDist)); 212 if (dmDist) { 213 CHKERRQ(DMDestroy(&dm)); 214 dm = dmDist; 215 } 216 CHKERRQ(PetscObjectSetName((PetscObject) dm, "Mesh")); 217 CHKERRQ(DMViewFromOptions(dm, NULL, "-dm_view")); 218 CHKERRQ(DMGetLabel(dm, name, &label)); 219 CHKERRQ(DMLabelView(label, PETSC_VIEWER_STDOUT_WORLD)); 220 CHKERRQ(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 CHKERRQ(PetscOptionsGetBool(NULL, NULL, "-universal", &run, NULL)); 234 if (!run) PetscFunctionReturn(0); 235 236 char filename[PETSC_MAX_PATH_LEN]; 237 PetscBool flg; 238 239 CHKERRQ(PetscOptionsGetString(NULL, NULL, "-filename", filename, sizeof(filename), &flg)); 240 if (flg) { 241 CHKERRQ(DMPlexCreateFromFile(comm, filename, "ex11_plex", PETSC_TRUE, &dm1)); 242 } else { 243 CHKERRQ(DMCreate(comm, &dm1)); 244 CHKERRQ(DMSetType(dm1, DMPLEX)); 245 CHKERRQ(DMSetFromOptions(dm1)); 246 } 247 CHKERRQ(DMHasLabel(dm1, "marker", ¬File)); 248 if (notFile) { 249 CHKERRQ(DMCreateLabel(dm1, "Boundary Faces")); 250 CHKERRQ(DMGetLabel(dm1, "Boundary Faces", &bd1)); 251 CHKERRQ(DMPlexMarkBoundaryFaces(dm1, 13, bd1)); 252 CHKERRQ(DMCreateLabel(dm1, "Boundary")); 253 CHKERRQ(DMGetLabel(dm1, "Boundary", &bd2)); 254 CHKERRQ(DMPlexMarkBoundaryFaces(dm1, 121, bd2)); 255 CHKERRQ(DMPlexLabelComplete(dm1, bd2)); 256 } 257 CHKERRQ(PetscObjectSetName((PetscObject) dm1, "First Mesh")); 258 CHKERRQ(DMViewFromOptions(dm1, NULL, "-dm_view")); 259 260 CHKERRQ(DMUniversalLabelCreate(dm1, &universal)); 261 CHKERRQ(DMUniversalLabelGetLabel(universal, &ulabel)); 262 CHKERRQ(PetscObjectViewFromOptions((PetscObject) ulabel, NULL, "-universal_view")); 263 264 if (!notFile) { 265 PetscInt Nl, l; 266 267 CHKERRQ(DMClone(dm1, &dm2)); 268 CHKERRQ(DMGetNumLabels(dm2, &Nl)); 269 for (l = Nl-1; l >= 0; --l) { 270 PetscBool isdepth, iscelltype; 271 const char *name; 272 273 CHKERRQ(DMGetLabelName(dm2, l, &name)); 274 CHKERRQ(PetscStrncmp(name, "depth", 6, &isdepth)); 275 CHKERRQ(PetscStrncmp(name, "celltype", 9, &iscelltype)); 276 if (!isdepth && !iscelltype) CHKERRQ(DMRemoveLabel(dm2, name, NULL)); 277 } 278 } else { 279 CHKERRQ(DMCreate(comm, &dm2)); 280 CHKERRQ(DMSetType(dm2, DMPLEX)); 281 CHKERRQ(DMSetFromOptions(dm2)); 282 } 283 CHKERRQ(PetscObjectSetName((PetscObject) dm2, "Second Mesh")); 284 CHKERRQ(DMUniversalLabelCreateLabels(universal, PETSC_TRUE, dm2)); 285 CHKERRQ(DMPlexGetChart(dm2, &pStart, &pEnd)); 286 for (p = pStart; p < pEnd; ++p) { 287 PetscInt val; 288 289 CHKERRQ(DMLabelGetValue(ulabel, p, &val)); 290 if (val < 0) continue; 291 CHKERRQ(DMUniversalLabelSetLabelValue(universal, dm2, PETSC_TRUE, p, val)); 292 } 293 CHKERRQ(DMViewFromOptions(dm2, NULL, "-dm_view")); 294 295 CHKERRQ(DMUniversalLabelDestroy(&universal)); 296 CHKERRQ(DMDestroy(&dm1)); 297 CHKERRQ(DMDestroy(&dm2)); 298 PetscFunctionReturn(0); 299 } 300 301 int main(int argc, char **argv) 302 { 303 304 CHKERRQ(PetscInitialize(&argc, &argv, NULL, help)); 305 /*CHKERRQ(ProcessOptions(PETSC_COMM_WORLD, &user));*/ 306 CHKERRQ(TestInsertion()); 307 CHKERRQ(TestEmptyStrata(PETSC_COMM_WORLD)); 308 CHKERRQ(TestDistribution(PETSC_COMM_WORLD)); 309 CHKERRQ(TestUniversalLabel(PETSC_COMM_WORLD)); 310 CHKERRQ(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