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