1c4762a1bSJed Brown static char help[] = "Create a Plex sphere from quads and create a P1 section\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscdmplex.h> 4c4762a1bSJed Brown 5c4762a1bSJed Brown static PetscErrorCode SetupSection(DM dm) 6c4762a1bSJed Brown { 7c4762a1bSJed Brown PetscSection s; 8c4762a1bSJed Brown PetscInt vStart, vEnd, v; 9c4762a1bSJed Brown 10c4762a1bSJed Brown PetscFunctionBeginUser; 119566063dSJacob Faibussowitsch PetscCall(DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd)); 129566063dSJacob Faibussowitsch PetscCall(PetscSectionCreate(PetscObjectComm((PetscObject) dm), &s)); 139566063dSJacob Faibussowitsch PetscCall(PetscSectionSetNumFields(s, 1)); 149566063dSJacob Faibussowitsch PetscCall(PetscSectionSetFieldComponents(s, 0, 1)); 159566063dSJacob Faibussowitsch PetscCall(PetscSectionSetChart(s, vStart, vEnd)); 16c4762a1bSJed Brown for (v = vStart; v < vEnd; ++v) { 179566063dSJacob Faibussowitsch PetscCall(PetscSectionSetDof(s, v, 1)); 189566063dSJacob Faibussowitsch PetscCall(PetscSectionSetFieldDof(s, v, 0, 1)); 19c4762a1bSJed Brown } 209566063dSJacob Faibussowitsch PetscCall(PetscSectionSetUp(s)); 219566063dSJacob Faibussowitsch PetscCall(DMSetLocalSection(dm, s)); 229566063dSJacob Faibussowitsch PetscCall(PetscSectionDestroy(&s)); 23c4762a1bSJed Brown PetscFunctionReturn(0); 24c4762a1bSJed Brown } 25c4762a1bSJed Brown 26c4762a1bSJed Brown int main(int argc, char **argv) 27c4762a1bSJed Brown { 28c4762a1bSJed Brown DM dm; 29c4762a1bSJed Brown Vec u; 30c4762a1bSJed Brown 31*327415f7SBarry Smith PetscFunctionBeginUser; 329566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &argv, NULL,help)); 339566063dSJacob Faibussowitsch PetscCall(DMCreate(PETSC_COMM_WORLD, &dm)); 349566063dSJacob Faibussowitsch PetscCall(DMSetType(dm, DMPLEX)); 359566063dSJacob Faibussowitsch PetscCall(DMSetFromOptions(dm)); 369566063dSJacob Faibussowitsch PetscCall(PetscObjectSetName((PetscObject) dm, "Sphere")); 379566063dSJacob Faibussowitsch PetscCall(DMViewFromOptions(dm, NULL, "-dm_view")); 3830602db0SMatthew G. Knepley 399566063dSJacob Faibussowitsch PetscCall(SetupSection(dm)); 409566063dSJacob Faibussowitsch PetscCall(DMGetGlobalVector(dm, &u)); 419566063dSJacob Faibussowitsch PetscCall(VecSet(u, 2)); 429566063dSJacob Faibussowitsch PetscCall(VecViewFromOptions(u, NULL, "-vec_view")); 439566063dSJacob Faibussowitsch PetscCall(DMRestoreGlobalVector(dm, &u)); 449566063dSJacob Faibussowitsch PetscCall(DMDestroy(&dm)); 459566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 46b122ec5aSJacob Faibussowitsch return 0; 47c4762a1bSJed Brown } 48c4762a1bSJed Brown 49c4762a1bSJed Brown /*TEST 50c4762a1bSJed Brown 5130602db0SMatthew G. Knepley testset: 52c4762a1bSJed Brown requires: !__float128 5330602db0SMatthew G. Knepley args: -dm_plex_shape sphere -dm_view 54c4762a1bSJed Brown 55c4762a1bSJed Brown test: 5630602db0SMatthew G. Knepley suffix: 2d_quad 5730602db0SMatthew G. Knepley args: -dm_plex_simplex 0 58c4762a1bSJed Brown 59c4762a1bSJed Brown test: 60c4762a1bSJed Brown suffix: 2d_tri 6130602db0SMatthew G. Knepley args: 62c4762a1bSJed Brown 63c4762a1bSJed Brown test: 64c4762a1bSJed Brown suffix: 3d_tri 6530602db0SMatthew G. Knepley args: -dm_plex_dim 3 6630602db0SMatthew G. Knepley 6730602db0SMatthew G. Knepley testset: 68c4762a1bSJed Brown requires: !__float128 69e600fa54SMatthew G. Knepley args: -dm_plex_shape sphere -petscpartitioner_type simple -dm_view 7030602db0SMatthew G. Knepley 7130602db0SMatthew G. Knepley test: 7230602db0SMatthew G. Knepley suffix: 2d_quad_parallel 7330602db0SMatthew G. Knepley nsize: 2 7430602db0SMatthew G. Knepley args: -dm_plex_simplex 0 7530602db0SMatthew G. Knepley 7630602db0SMatthew G. Knepley test: 7730602db0SMatthew G. Knepley suffix: 2d_tri_parallel 7830602db0SMatthew G. Knepley nsize: 2 79c4762a1bSJed Brown 80c4762a1bSJed Brown test: 81c4762a1bSJed Brown suffix: 3d_tri_parallel 82c4762a1bSJed Brown nsize: 2 8330602db0SMatthew G. Knepley args: -dm_plex_dim 3 84c4762a1bSJed Brown 85c4762a1bSJed Brown TEST*/ 86