1c4762a1bSJed Brown static char help[] = "Define a simple field over the mesh\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscdmplex.h> 4c4762a1bSJed Brown 5c4762a1bSJed Brown int main(int argc, char **argv) 6c4762a1bSJed Brown { 7*30602db0SMatthew G. Knepley DM dm; 8c4762a1bSJed Brown Vec u; 9c4762a1bSJed Brown PetscSection section; 10c4762a1bSJed Brown PetscViewer viewer; 11*30602db0SMatthew G. Knepley PetscInt dim, numFields, numBC, i; 12c4762a1bSJed Brown PetscInt numComp[3]; 13c4762a1bSJed Brown PetscInt numDof[12]; 14c4762a1bSJed Brown PetscInt bcField[1]; 15c4762a1bSJed Brown IS bcPointIS[1]; 16c4762a1bSJed Brown PetscErrorCode ierr; 17c4762a1bSJed Brown 18c4762a1bSJed Brown ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 19c4762a1bSJed Brown /* Create a mesh */ 20*30602db0SMatthew G. Knepley ierr = DMCreate(PETSC_COMM_WORLD, &dm);CHKERRQ(ierr); 21*30602db0SMatthew G. Knepley ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 22*30602db0SMatthew G. Knepley ierr = DMSetFromOptions(dm);CHKERRQ(ierr); 23*30602db0SMatthew G. Knepley ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 24*30602db0SMatthew G. Knepley ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 25c4762a1bSJed Brown /* Create a scalar field u, a vector field v, and a surface vector field w */ 26c4762a1bSJed Brown numFields = 3; 27c4762a1bSJed Brown numComp[0] = 1; 28c4762a1bSJed Brown numComp[1] = dim; 29c4762a1bSJed Brown numComp[2] = dim-1; 30c4762a1bSJed Brown for (i = 0; i < numFields*(dim+1); ++i) numDof[i] = 0; 31c4762a1bSJed Brown /* Let u be defined on vertices */ 32c4762a1bSJed Brown numDof[0*(dim+1)+0] = 1; 33c4762a1bSJed Brown /* Let v be defined on cells */ 34c4762a1bSJed Brown numDof[1*(dim+1)+dim] = dim; 35c4762a1bSJed Brown /* Let w be defined on faces */ 36c4762a1bSJed Brown numDof[2*(dim+1)+dim-1] = dim-1; 37c4762a1bSJed Brown /* Setup boundary conditions */ 38c4762a1bSJed Brown numBC = 1; 39c4762a1bSJed Brown /* Prescribe a Dirichlet condition on u on the boundary 40c4762a1bSJed Brown Label "marker" is made by the mesh creation routine */ 41c4762a1bSJed Brown bcField[0] = 0; 42c4762a1bSJed Brown ierr = DMGetStratumIS(dm, "marker", 1, &bcPointIS[0]);CHKERRQ(ierr); 43c4762a1bSJed Brown /* Create a PetscSection with this data layout */ 44c4762a1bSJed Brown ierr = DMSetNumFields(dm, numFields);CHKERRQ(ierr); 45c4762a1bSJed Brown ierr = DMPlexCreateSection(dm, NULL, numComp, numDof, numBC, bcField, NULL, bcPointIS, NULL, §ion);CHKERRQ(ierr); 46c4762a1bSJed Brown ierr = ISDestroy(&bcPointIS[0]);CHKERRQ(ierr); 47c4762a1bSJed Brown /* Name the Field variables */ 48c4762a1bSJed Brown ierr = PetscSectionSetFieldName(section, 0, "u");CHKERRQ(ierr); 49c4762a1bSJed Brown ierr = PetscSectionSetFieldName(section, 1, "v");CHKERRQ(ierr); 50c4762a1bSJed Brown ierr = PetscSectionSetFieldName(section, 2, "w");CHKERRQ(ierr); 51c4762a1bSJed Brown ierr = PetscSectionView(section, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 52c4762a1bSJed Brown /* Tell the DM to use this data layout */ 53c4762a1bSJed Brown ierr = DMSetLocalSection(dm, section);CHKERRQ(ierr); 54c4762a1bSJed Brown /* Create a Vec with this layout and view it */ 55c4762a1bSJed Brown ierr = DMGetGlobalVector(dm, &u);CHKERRQ(ierr); 56c4762a1bSJed Brown ierr = PetscViewerCreate(PETSC_COMM_WORLD, &viewer);CHKERRQ(ierr); 57c4762a1bSJed Brown ierr = PetscViewerSetType(viewer, PETSCVIEWERVTK);CHKERRQ(ierr); 588ec8862eSJed Brown ierr = PetscViewerFileSetName(viewer, "sol.vtu");CHKERRQ(ierr); 59c4762a1bSJed Brown ierr = VecView(u, viewer);CHKERRQ(ierr); 60c4762a1bSJed Brown ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 61c4762a1bSJed Brown ierr = DMRestoreGlobalVector(dm, &u);CHKERRQ(ierr); 62c4762a1bSJed Brown /* Cleanup */ 63c4762a1bSJed Brown ierr = PetscSectionDestroy(§ion);CHKERRQ(ierr); 64c4762a1bSJed Brown ierr = DMDestroy(&dm);CHKERRQ(ierr); 65c4762a1bSJed Brown ierr = PetscFinalize(); 66c4762a1bSJed Brown return ierr; 67c4762a1bSJed Brown } 68c4762a1bSJed Brown 69c4762a1bSJed Brown /*TEST 70c4762a1bSJed Brown 71c4762a1bSJed Brown test: 72c4762a1bSJed Brown suffix: 0 73c4762a1bSJed Brown requires: triangle 74*30602db0SMatthew G. Knepley args: -info :~sys,mat 75c4762a1bSJed Brown test: 76c4762a1bSJed Brown suffix: 1 77c4762a1bSJed Brown requires: ctetgen 78*30602db0SMatthew G. Knepley args: -dm_plex_dim 3 -info :~sys,mat 79c4762a1bSJed Brown 80c4762a1bSJed Brown TEST*/ 81