1c4762a1bSJed Brown static char help[] = "Tests for point location\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscsf.h> 4c4762a1bSJed Brown #include <petscdmplex.h> 5c4762a1bSJed Brown 630602db0SMatthew G. Knepley static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm) 7c4762a1bSJed Brown { 8c4762a1bSJed Brown PetscErrorCode ierr; 9c4762a1bSJed Brown 10c4762a1bSJed Brown PetscFunctionBeginUser; 1130602db0SMatthew G. Knepley ierr = DMCreate(comm, dm);CHKERRQ(ierr); 1230602db0SMatthew G. Knepley ierr = DMSetType(*dm, DMPLEX);CHKERRQ(ierr); 1330602db0SMatthew G. Knepley ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 14c4762a1bSJed Brown ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr); 15c4762a1bSJed Brown PetscFunctionReturn(0); 16c4762a1bSJed Brown } 17c4762a1bSJed Brown 1830602db0SMatthew G. Knepley static PetscErrorCode TestLocation(DM dm) 19c4762a1bSJed Brown { 2030602db0SMatthew G. Knepley PetscInt dim; 21c4762a1bSJed Brown PetscInt cStart, cEnd, c; 22c4762a1bSJed Brown PetscErrorCode ierr; 23c4762a1bSJed Brown 24c4762a1bSJed Brown PetscFunctionBeginUser; 2530602db0SMatthew G. Knepley ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 26c4762a1bSJed Brown ierr = DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);CHKERRQ(ierr); 27c4762a1bSJed Brown /* Locate all centroids */ 28c4762a1bSJed Brown for (c = cStart; c < cEnd; ++c) { 29c4762a1bSJed Brown Vec v; 30c4762a1bSJed Brown PetscSF cellSF = NULL; 31c4762a1bSJed Brown const PetscSFNode *cells; 32c4762a1bSJed Brown PetscScalar *a; 33c4762a1bSJed Brown PetscReal centroid[3]; 34c4762a1bSJed Brown PetscInt n, d; 35c4762a1bSJed Brown 36c4762a1bSJed Brown ierr = DMPlexComputeCellGeometryFVM(dm, c, NULL, centroid, NULL);CHKERRQ(ierr); 37c4762a1bSJed Brown ierr = VecCreateSeq(PETSC_COMM_SELF, dim, &v);CHKERRQ(ierr); 38c4762a1bSJed Brown ierr = VecSetBlockSize(v, dim);CHKERRQ(ierr); 39c4762a1bSJed Brown ierr = VecGetArray(v, &a);CHKERRQ(ierr); 40c4762a1bSJed Brown for (d = 0; d < dim; ++d) a[d] = centroid[d]; 41c4762a1bSJed Brown ierr = VecRestoreArray(v, &a);CHKERRQ(ierr); 42c4762a1bSJed Brown ierr = DMLocatePoints(dm, v, DM_POINTLOCATION_NONE, &cellSF);CHKERRQ(ierr); 43c4762a1bSJed Brown ierr = VecDestroy(&v);CHKERRQ(ierr); 44c4762a1bSJed Brown ierr = PetscSFGetGraph(cellSF,NULL,&n,NULL,&cells);CHKERRQ(ierr); 45*1af33867SMatthew G. Knepley if (n != 1) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Cell %D: Found %d cells instead of 1", c, n); 46*1af33867SMatthew G. Knepley if (cells[0].index != c) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Could not locate centroid of cell %D, instead found %D", c, cells[0].index); 47c4762a1bSJed Brown ierr = PetscSFDestroy(&cellSF);CHKERRQ(ierr); 48c4762a1bSJed Brown } 49c4762a1bSJed Brown PetscFunctionReturn(0); 50c4762a1bSJed Brown } 51c4762a1bSJed Brown 52c4762a1bSJed Brown int main(int argc, char **argv) 53c4762a1bSJed Brown { 54c4762a1bSJed Brown DM dm; 55c4762a1bSJed Brown PetscErrorCode ierr; 56c4762a1bSJed Brown 57c4762a1bSJed Brown ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 5830602db0SMatthew G. Knepley ierr = CreateMesh(PETSC_COMM_WORLD, &dm);CHKERRQ(ierr); 5930602db0SMatthew G. Knepley ierr = TestLocation(dm);CHKERRQ(ierr); 60c4762a1bSJed Brown ierr = DMDestroy(&dm);CHKERRQ(ierr); 61c4762a1bSJed Brown ierr = PetscFinalize(); 62c4762a1bSJed Brown return ierr; 63c4762a1bSJed Brown } 64c4762a1bSJed Brown 65c4762a1bSJed Brown /*TEST 66c4762a1bSJed Brown 67c4762a1bSJed Brown test: 68*1af33867SMatthew G. Knepley suffix: seg 69*1af33867SMatthew G. Knepley args: -dm_plex_dim 1 -dm_plex_box_faces 10 70*1af33867SMatthew G. Knepley 71*1af33867SMatthew G. Knepley test: 72*1af33867SMatthew G. Knepley suffix: tri 73c4762a1bSJed Brown requires: triangle 74*1af33867SMatthew G. Knepley args: -dm_coord_space 0 -dm_plex_box_faces 5,5 75*1af33867SMatthew G. Knepley 76*1af33867SMatthew G. Knepley test: 77*1af33867SMatthew G. Knepley suffix: quad 78*1af33867SMatthew G. Knepley args: -dm_plex_simplex 0 -dm_plex_box_faces 5,5 79*1af33867SMatthew G. Knepley 80*1af33867SMatthew G. Knepley test: 81*1af33867SMatthew G. Knepley suffix: tet 82*1af33867SMatthew G. Knepley requires: ctetgen 83*1af33867SMatthew G. Knepley args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 84*1af33867SMatthew G. Knepley 85*1af33867SMatthew G. Knepley test: 86*1af33867SMatthew G. Knepley suffix: hex 87*1af33867SMatthew G. Knepley args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 3,3,3 88c4762a1bSJed Brown 89c4762a1bSJed Brown TEST*/ 90