1c4762a1bSJed Brown program DMPlexTestField 2c4762a1bSJed Brown#include "petsc/finclude/petscdmplex.h" 3c4762a1bSJed Brown#include "petsc/finclude/petscdmlabel.h" 4c4762a1bSJed Brown use petscdmplex 5c4762a1bSJed Brown use petscsys 6c4762a1bSJed Brown implicit none 7c4762a1bSJed Brown 8c4762a1bSJed Brown DM :: dm 9c4762a1bSJed Brown DMLabel :: label 10c4762a1bSJed Brown Vec :: u 11c4762a1bSJed Brown PetscViewer :: viewer 12c4762a1bSJed Brown PetscSection :: section 13c4762a1bSJed Brown PetscInt :: dim,numFields,numBC 14c4762a1bSJed Brown PetscInt :: i,val 15c4762a1bSJed Brown DMLabel, pointer :: nolabel(:) => NULL() 16c4762a1bSJed Brown PetscInt, target, dimension(3) :: numComp 17c4762a1bSJed Brown PetscInt, pointer :: pNumComp(:) 18c4762a1bSJed Brown PetscInt, target, dimension(12) :: numDof 19c4762a1bSJed Brown PetscInt, pointer :: pNumDof(:) 20c4762a1bSJed Brown PetscInt, target, dimension(1) :: bcField 21c4762a1bSJed Brown PetscInt, pointer :: pBcField(:) 22c4762a1bSJed Brown PetscInt, parameter :: zero = 0, one = 1, two = 2, eight = 8 232e3d3ef9SMatthew G. Knepley PetscMPIInt :: size 24c4762a1bSJed Brown IS, target, dimension(1) :: bcCompIS 25c4762a1bSJed Brown IS, target, dimension(1) :: bcPointIS 26c4762a1bSJed Brown IS, pointer :: pBcCompIS(:) 27c4762a1bSJed Brown IS, pointer :: pBcPointIS(:) 28c4762a1bSJed Brown PetscErrorCode :: ierr 29c4762a1bSJed Brown 30*d8606c27SBarry Smith PetscCallA(PetscInitialize(ierr)) 31*d8606c27SBarry Smith PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD, size, ierr)) 32c4762a1bSJed Brown! Create a mesh 33*d8606c27SBarry Smith PetscCallA(DMCreate(PETSC_COMM_WORLD, dm, ierr)) 34*d8606c27SBarry Smith PetscCallA(DMSetType(dm, DMPLEX, ierr)) 35*d8606c27SBarry Smith PetscCallA(DMSetFromOptions(dm, ierr)) 36*d8606c27SBarry Smith PetscCallA(DMViewFromOptions(dm, PETSC_NULL_VEC, '-dm_view', ierr)) 37*d8606c27SBarry Smith PetscCallA(DMGetDimension(dm, dim, ierr)) 38c4762a1bSJed Brown! Create a scalar field u, a vector field v, and a surface vector field w 39c4762a1bSJed Brown numFields = 3 40c4762a1bSJed Brown numComp(1) = 1 41c4762a1bSJed Brown numComp(2) = dim 42c4762a1bSJed Brown numComp(3) = dim-1 43c4762a1bSJed Brown pNumComp => numComp 44c4762a1bSJed Brown do i = 1, numFields*(dim+1) 45c4762a1bSJed Brown numDof(i) = 0 46c4762a1bSJed Brown end do 47c4762a1bSJed Brown! Let u be defined on vertices 48c4762a1bSJed Brown numDof(0*(dim+1)+1) = 1 49c4762a1bSJed Brown! Let v be defined on cells 50c4762a1bSJed Brown numDof(1*(dim+1)+dim+1) = dim 51c4762a1bSJed Brown! Let v be defined on faces 52c4762a1bSJed Brown numDof(2*(dim+1)+dim) = dim-1 53c4762a1bSJed Brown pNumDof => numDof 54c4762a1bSJed Brown! Setup boundary conditions 55c4762a1bSJed Brown numBC = 1 56c4762a1bSJed Brown! Test label retrieval 57*d8606c27SBarry Smith PetscCallA(DMGetLabel(dm, 'marker', label, ierr)) 58*d8606c27SBarry Smith PetscCallA(DMLabelGetValue(label, zero, val, ierr)) 592e3d3ef9SMatthew G. Knepley if (size .eq. 1 .and. val .ne. -1) then 60c4762a1bSJed Brown SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error in library') 61c4762a1bSJed Brown endif 62*d8606c27SBarry Smith PetscCallA(DMLabelGetValue(label, eight, val, ierr)) 632e3d3ef9SMatthew G. Knepley if (size .eq. 1 .and. val .ne. 1) then 64c4762a1bSJed Brown SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error in library') 65c4762a1bSJed Brown endif 66c4762a1bSJed Brown! Prescribe a Dirichlet condition on u on the boundary 67c4762a1bSJed Brown! Label "marker" is made by the mesh creation routine 68c4762a1bSJed Brown bcField(1) = 0 69c4762a1bSJed Brown pBcField => bcField 70*d8606c27SBarry Smith PetscCallA(ISCreateStride(PETSC_COMM_WORLD, one, zero, one, bcCompIS(1), ierr)) 71c4762a1bSJed Brown pBcCompIS => bcCompIS 72*d8606c27SBarry Smith PetscCallA(DMGetStratumIS(dm, 'marker', one, bcPointIS(1),ierr)) 73c4762a1bSJed Brown pBcPointIS => bcPointIS 74c4762a1bSJed Brown! Create a PetscSection with this data layout 75*d8606c27SBarry Smith PetscCallA(DMSetNumFields(dm, numFields,ierr)) 76*d8606c27SBarry Smith PetscCallA(DMPlexCreateSection(dm,nolabel,pNumComp,pNumDof,numBC,pBcField,pBcCompIS,pBcPointIS,PETSC_NULL_IS,section,ierr)) 77*d8606c27SBarry Smith PetscCallA(ISDestroy(bcCompIS(1), ierr)) 78*d8606c27SBarry Smith PetscCallA(ISDestroy(bcPointIS(1), ierr)) 79c4762a1bSJed Brown! Name the Field variables 80*d8606c27SBarry Smith PetscCallA(PetscSectionSetFieldName(section, zero, 'u', ierr)) 81*d8606c27SBarry Smith PetscCallA(PetscSectionSetFieldName(section, one, 'v', ierr)) 82*d8606c27SBarry Smith PetscCallA(PetscSectionSetFieldName(section, two, 'w', ierr)) 832e3d3ef9SMatthew G. Knepley if (size .eq. 1) then 84*d8606c27SBarry Smith PetscCallA(PetscSectionView(section, PETSC_VIEWER_STDOUT_WORLD, ierr)) 852e3d3ef9SMatthew G. Knepley endif 86c4762a1bSJed Brown! Tell the DM to use this data layout 87*d8606c27SBarry Smith PetscCallA(DMSetLocalSection(dm, section, ierr)) 88c4762a1bSJed Brown! Create a Vec with this layout and view it 89*d8606c27SBarry Smith PetscCallA(DMGetGlobalVector(dm, u, ierr)) 90*d8606c27SBarry Smith PetscCallA(PetscViewerCreate(PETSC_COMM_WORLD, viewer, ierr)) 91*d8606c27SBarry Smith PetscCallA(PetscViewerSetType(viewer, PETSCVIEWERVTK, ierr)) 92*d8606c27SBarry Smith PetscCallA(PetscViewerFileSetName(viewer, 'sol.vtu', ierr)) 93*d8606c27SBarry Smith PetscCallA(VecView(u, viewer, ierr)) 94*d8606c27SBarry Smith PetscCallA(PetscViewerDestroy(viewer, ierr)) 95*d8606c27SBarry Smith PetscCallA(DMRestoreGlobalVector(dm, u, ierr)) 96c4762a1bSJed Brown! Cleanup 97*d8606c27SBarry Smith PetscCallA(PetscSectionDestroy(section, ierr)) 98*d8606c27SBarry Smith PetscCallA(DMDestroy(dm, ierr)) 99c4762a1bSJed Brown 100*d8606c27SBarry Smith PetscCallA(PetscFinalize(ierr)) 101c4762a1bSJed Brown end program DMPlexTestField 102c4762a1bSJed Brown 103c4762a1bSJed Brown!/*TEST 104c4762a1bSJed Brown! build: 105dfd57a17SPierre Jolivet! requires: defined(PETSC_USING_F90FREEFORM) 106c4762a1bSJed Brown! 107c4762a1bSJed Brown! test: 108c4762a1bSJed Brown! suffix: 0 109c4762a1bSJed Brown! requires: triangle 11030602db0SMatthew G. Knepley! args: -info :~sys,mat: 111c4762a1bSJed Brown! 112c4762a1bSJed Brown! test: 1132e3d3ef9SMatthew G. Knepley! suffix: 0_2 1142e3d3ef9SMatthew G. Knepley! requires: triangle 1152e3d3ef9SMatthew G. Knepley! nsize: 2 1162e3d3ef9SMatthew G. Knepley! args: -info :~sys,mat,partitioner: 1172e3d3ef9SMatthew G. Knepley! 1182e3d3ef9SMatthew G. Knepley! test: 119c4762a1bSJed Brown! suffix: 1 120c4762a1bSJed Brown! requires: ctetgen 12130602db0SMatthew G. Knepley! args: -dm_plex_dim 3 -info :~sys,mat: 122c4762a1bSJed Brown! 123c4762a1bSJed Brown!TEST*/ 124