Actual source code: ex17.c
1: static char help[] = "Tests for point location\n\n";
3: #include <petscsf.h>
4: #include <petscdmplex.h>
6: static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm)
7: {
9: DMCreate(comm, dm);
10: DMSetType(*dm, DMPLEX);
11: DMSetFromOptions(*dm);
12: DMViewFromOptions(*dm, NULL, "-dm_view");
13: return 0;
14: }
16: static PetscErrorCode TestLocation(DM dm)
17: {
18: Vec points;
19: PetscSF cellSF = NULL;
20: const PetscSFNode *cells;
21: PetscScalar *a;
22: PetscInt cdim, n;
23: PetscInt cStart, cEnd, c;
26: DMGetCoordinateDim(dm, &cdim);
27: DMGetCoordinatesLocalSetUp(dm);
28: DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);
29: /* Locate all centroids */
30: VecCreateSeq(PETSC_COMM_SELF, (cEnd - cStart) * cdim, &points);
31: VecSetBlockSize(points, cdim);
32: VecGetArray(points, &a);
33: for (c = cStart; c < cEnd; ++c) {
34: PetscReal centroid[3];
35: PetscInt off = (c - cStart) * cdim, d;
37: DMPlexComputeCellGeometryFVM(dm, c, NULL, centroid, NULL);
38: for (d = 0; d < cdim; ++d) a[off + d] = centroid[d];
39: }
40: VecRestoreArray(points, &a);
41: DMLocatePoints(dm, points, DM_POINTLOCATION_NONE, &cellSF);
42: VecDestroy(&points);
43: PetscSFGetGraph(cellSF, NULL, &n, NULL, &cells);
44: if (n != (cEnd - cStart)) {
45: for (c = 0; c < n; ++c) {
46: if (cells[c].index != c + cStart) PetscPrintf(PETSC_COMM_SELF, "Could not locate centroid of cell %" PetscInt_FMT ", error %" PetscInt_FMT "\n", c + cStart, cells[c].index);
47: }
48: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Located %" PetscInt_FMT " points instead of %" PetscInt_FMT, n, cEnd - cStart);
49: }
51: PetscSFDestroy(&cellSF);
52: return 0;
53: }
55: int main(int argc, char **argv)
56: {
57: DM dm;
60: PetscInitialize(&argc, &argv, NULL, help);
61: CreateMesh(PETSC_COMM_WORLD, &dm);
62: TestLocation(dm);
63: DMDestroy(&dm);
64: PetscFinalize();
65: return 0;
66: }
68: /*TEST
70: testset:
71: args: -dm_plex_dim 1 -dm_plex_box_faces 10
73: test:
74: suffix: seg
76: test:
77: suffix: seg_hash
78: args: -dm_refine 2 -dm_plex_hash_location
80: testset:
81: args: -dm_plex_box_faces 5,5
83: test:
84: suffix: tri
85: requires: triangle
87: test:
88: suffix: tri_hash
89: requires: triangle
90: args: -dm_refine 2 -dm_plex_hash_location
92: test:
93: suffix: quad
94: args: -dm_plex_simplex 0
96: test:
97: suffix: quad_hash
98: args: -dm_plex_simplex 0 -dm_refine 2 -dm_plex_hash_location
100: testset:
101: args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3
103: test:
104: suffix: tet
105: requires: ctetgen
107: test:
108: suffix: tet_hash
109: requires: ctetgen
110: args: -dm_refine 1 -dm_plex_hash_location
112: test:
113: suffix: hex
114: args: -dm_plex_simplex 0
116: test:
117: suffix: hex_hash
118: args: -dm_plex_simplex 0 -dm_refine 1 -dm_plex_hash_location
120: TEST*/