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 23c4762a1bSJed Brown IS, target, dimension(1) :: bcCompIS 24c4762a1bSJed Brown IS, target, dimension(1) :: bcPointIS 25c4762a1bSJed Brown IS, pointer :: pBcCompIS(:) 26c4762a1bSJed Brown IS, pointer :: pBcPointIS(:) 27c4762a1bSJed Brown PetscErrorCode :: ierr 28c4762a1bSJed Brown 29c4762a1bSJed Brown call PetscInitialize(PETSC_NULL_CHARACTER, ierr) 30c4762a1bSJed Brown if (ierr .ne. 0) then 31c4762a1bSJed Brown print*,'Unable to initialize PETSc' 32c4762a1bSJed Brown stop 33c4762a1bSJed Brown endif 34c4762a1bSJed Brown! Create a mesh 3530602db0SMatthew G. Knepley call DMCreate(PETSC_COMM_WORLD, dm, ierr);CHKERRA(ierr) 3630602db0SMatthew G. Knepley call DMSetType(dm, DMPLEX, ierr);CHKERRA(ierr) 3730602db0SMatthew G. Knepley call DMSetFromOptions(dm, ierr);CHKERRA(ierr) 3830602db0SMatthew G. Knepley call DMGetDimension(dm, dim, ierr);CHKERRA(ierr) 39c4762a1bSJed Brown! Create a scalar field u, a vector field v, and a surface vector field w 40c4762a1bSJed Brown numFields = 3 41c4762a1bSJed Brown numComp(1) = 1 42c4762a1bSJed Brown numComp(2) = dim 43c4762a1bSJed Brown numComp(3) = dim-1 44c4762a1bSJed Brown pNumComp => numComp 45c4762a1bSJed Brown do i = 1, numFields*(dim+1) 46c4762a1bSJed Brown numDof(i) = 0 47c4762a1bSJed Brown end do 48c4762a1bSJed Brown! Let u be defined on vertices 49c4762a1bSJed Brown numDof(0*(dim+1)+1) = 1 50c4762a1bSJed Brown! Let v be defined on cells 51c4762a1bSJed Brown numDof(1*(dim+1)+dim+1) = dim 52c4762a1bSJed Brown! Let v be defined on faces 53c4762a1bSJed Brown numDof(2*(dim+1)+dim) = dim-1 54c4762a1bSJed Brown pNumDof => numDof 55c4762a1bSJed Brown! Setup boundary conditions 56c4762a1bSJed Brown numBC = 1 57c4762a1bSJed Brown! Test label retrieval 58c4762a1bSJed Brown call DMGetLabel(dm, 'marker', label, ierr);CHKERRA(ierr) 59c4762a1bSJed Brown call DMLabelGetValue(label, zero, val, ierr);CHKERRA(ierr) 60c4762a1bSJed Brown if (val .ne. -1) then 61c4762a1bSJed Brown SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error in library') 62c4762a1bSJed Brown endif 63c4762a1bSJed Brown call DMLabelGetValue(label, eight, val, ierr);CHKERRA(ierr) 64c4762a1bSJed Brown if (val .ne. 1) then 65c4762a1bSJed Brown SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error in library') 66c4762a1bSJed Brown endif 67c4762a1bSJed Brown! Prescribe a Dirichlet condition on u on the boundary 68c4762a1bSJed Brown! Label "marker" is made by the mesh creation routine 69c4762a1bSJed Brown bcField(1) = 0 70c4762a1bSJed Brown pBcField => bcField 71c4762a1bSJed Brown call ISCreateStride(PETSC_COMM_WORLD, one, zero, one, bcCompIS(1), ierr);CHKERRA(ierr) 72c4762a1bSJed Brown pBcCompIS => bcCompIS 73c4762a1bSJed Brown call DMGetStratumIS(dm, 'marker', one, bcPointIS(1),ierr);CHKERRA(ierr) 74c4762a1bSJed Brown pBcPointIS => bcPointIS 75c4762a1bSJed Brown! Create a PetscSection with this data layout 76c4762a1bSJed Brown call DMSetNumFields(dm, numFields,ierr);CHKERRA(ierr) 77c4762a1bSJed Brown call DMPlexCreateSection(dm,nolabel,pNumComp,pNumDof,numBC,pBcField,pBcCompIS,pBcPointIS,PETSC_NULL_IS,section,ierr) 78c4762a1bSJed Brown CHKERRA(ierr) 79c4762a1bSJed Brown call ISDestroy(bcCompIS(1), ierr);CHKERRA(ierr) 80c4762a1bSJed Brown call ISDestroy(bcPointIS(1), ierr);CHKERRA(ierr) 81c4762a1bSJed Brown! Name the Field variables 82c4762a1bSJed Brown call PetscSectionSetFieldName(section, zero, 'u', ierr);CHKERRA(ierr) 83c4762a1bSJed Brown call PetscSectionSetFieldName(section, one, 'v', ierr);CHKERRA(ierr) 84c4762a1bSJed Brown call PetscSectionSetFieldName(section, two, 'w', ierr);CHKERRA(ierr) 85c4762a1bSJed Brown call PetscSectionView(section, PETSC_VIEWER_STDOUT_WORLD, ierr);CHKERRA(ierr) 86c4762a1bSJed Brown! Tell the DM to use this data layout 87c4762a1bSJed Brown call DMSetLocalSection(dm, section, ierr);CHKERRA(ierr) 88c4762a1bSJed Brown! Create a Vec with this layout and view it 89c4762a1bSJed Brown call DMGetGlobalVector(dm, u, ierr);CHKERRA(ierr) 90c4762a1bSJed Brown call PetscViewerCreate(PETSC_COMM_WORLD, viewer, ierr);CHKERRA(ierr) 91c4762a1bSJed Brown call PetscViewerSetType(viewer, PETSCVIEWERVTK, ierr);CHKERRA(ierr) 928ec8862eSJed Brown call PetscViewerFileSetName(viewer, 'sol.vtu', ierr);CHKERRA(ierr) 93c4762a1bSJed Brown call VecView(u, viewer, ierr);CHKERRA(ierr) 94c4762a1bSJed Brown call PetscViewerDestroy(viewer, ierr);CHKERRA(ierr) 95c4762a1bSJed Brown call DMRestoreGlobalVector(dm, u, ierr);CHKERRA(ierr) 96c4762a1bSJed Brown! Cleanup 97c4762a1bSJed Brown call PetscSectionDestroy(section, ierr);CHKERRA(ierr) 98c4762a1bSJed Brown call DMDestroy(dm, ierr);CHKERRA(ierr) 99c4762a1bSJed Brown 100c4762a1bSJed Brown call PetscFinalize(ierr) 101c4762a1bSJed Brown end program DMPlexTestField 102c4762a1bSJed Brown 103c4762a1bSJed Brown!/*TEST 104c4762a1bSJed Brown! build: 105*dfd57a17SPierre 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: 113c4762a1bSJed Brown! suffix: 1 114c4762a1bSJed Brown! requires: ctetgen 11530602db0SMatthew G. Knepley! args: -dm_plex_dim 3 -info :~sys,mat: 116c4762a1bSJed Brown! 117c4762a1bSJed Brown!TEST*/ 118