1ce78bad3SBarry Smith#include <petsc/finclude/petscdmplex.h> 2ce78bad3SBarry Smith#include <petsc/finclude/petscdmlabel.h> 3c5e229c2SMartin Diehlprogram DMPlexTestField 4ce78bad3SBarry Smith use petscdm 501fa2b5aSMartin Diehl use petscdmplex 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 PetscInt, target, dimension(3) :: numComp 16c4762a1bSJed Brown PetscInt, pointer :: pNumComp(:) 17c4762a1bSJed Brown PetscInt, target, dimension(12) :: numDof 18c4762a1bSJed Brown PetscInt, pointer :: pNumDof(:) 19c4762a1bSJed Brown PetscInt, target, dimension(1) :: bcField 20c4762a1bSJed Brown PetscInt, pointer :: pBcField(:) 21c4762a1bSJed Brown PetscInt, parameter :: zero = 0, one = 1, two = 2, eight = 8 222e3d3ef9SMatthew G. Knepley PetscMPIInt :: size 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 29d8606c27SBarry Smith PetscCallA(PetscInitialize(ierr)) 30d8606c27SBarry Smith PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD, size, ierr)) 31c4762a1bSJed Brown! Create a mesh 32d8606c27SBarry Smith PetscCallA(DMCreate(PETSC_COMM_WORLD, dm, ierr)) 33d8606c27SBarry Smith PetscCallA(DMSetType(dm, DMPLEX, ierr)) 34d8606c27SBarry Smith PetscCallA(DMSetFromOptions(dm, ierr)) 35ce78bad3SBarry Smith PetscCallA(DMViewFromOptions(dm, PETSC_NULL_OBJECT, '-dm_view', ierr)) 36d8606c27SBarry Smith PetscCallA(DMGetDimension(dm, dim, ierr)) 37c4762a1bSJed Brown! Create a scalar field u, a vector field v, and a surface vector field w 38c4762a1bSJed Brown numFields = 3 39c4762a1bSJed Brown numComp(1) = 1 40c4762a1bSJed Brown numComp(2) = dim 41c4762a1bSJed Brown numComp(3) = dim - 1 42c4762a1bSJed Brown pNumComp => numComp 43c4762a1bSJed Brown do i = 1, numFields*(dim + 1) 44c4762a1bSJed Brown numDof(i) = 0 45c4762a1bSJed Brown end do 46c4762a1bSJed Brown! Let u be defined on vertices 47c4762a1bSJed Brown numDof(0*(dim + 1) + 1) = 1 48c4762a1bSJed Brown! Let v be defined on cells 49c4762a1bSJed Brown numDof(1*(dim + 1) + dim + 1) = dim 50c4762a1bSJed Brown! Let v be defined on faces 51c4762a1bSJed Brown numDof(2*(dim + 1) + dim) = dim - 1 52c4762a1bSJed Brown pNumDof => numDof 53c4762a1bSJed Brown! Setup boundary conditions 54c4762a1bSJed Brown numBC = 1 55c4762a1bSJed Brown! Test label retrieval 56d8606c27SBarry Smith PetscCallA(DMGetLabel(dm, 'marker', label, ierr)) 57d8606c27SBarry Smith PetscCallA(DMLabelGetValue(label, zero, val, ierr)) 584820e4eaSBarry Smith PetscCheckA(size /= 1 .or. val == -1, PETSC_COMM_SELF, PETSC_ERR_PLIB, 'Error in library') 59d8606c27SBarry Smith PetscCallA(DMLabelGetValue(label, eight, val, ierr)) 604820e4eaSBarry Smith PetscCheckA(size /= 1 .or. val == 1, PETSC_COMM_SELF, PETSC_ERR_PLIB, 'Error in library') 61c4762a1bSJed Brown! Prescribe a Dirichlet condition on u on the boundary 62c4762a1bSJed Brown! Label "marker" is made by the mesh creation routine 63c4762a1bSJed Brown bcField(1) = 0 64c4762a1bSJed Brown pBcField => bcField 65d8606c27SBarry Smith PetscCallA(ISCreateStride(PETSC_COMM_WORLD, one, zero, one, bcCompIS(1), ierr)) 66c4762a1bSJed Brown pBcCompIS => bcCompIS 67d8606c27SBarry Smith PetscCallA(DMGetStratumIS(dm, 'marker', one, bcPointIS(1), ierr)) 68c4762a1bSJed Brown pBcPointIS => bcPointIS 69c4762a1bSJed Brown! Create a PetscSection with this data layout 70d8606c27SBarry Smith PetscCallA(DMSetNumFields(dm, numFields, ierr)) 71ce78bad3SBarry Smith PetscCallA(DMPlexCreateSection(dm, PETSC_NULL_DMLABEL_ARRAY, pNumComp, pNumDof, numBC, pBcField, pBcCompIS, pBcPointIS, PETSC_NULL_IS, section, ierr)) 72d8606c27SBarry Smith PetscCallA(ISDestroy(bcCompIS(1), ierr)) 73d8606c27SBarry Smith PetscCallA(ISDestroy(bcPointIS(1), ierr)) 74c4762a1bSJed Brown! Name the Field variables 75d8606c27SBarry Smith PetscCallA(PetscSectionSetFieldName(section, zero, 'u', ierr)) 76d8606c27SBarry Smith PetscCallA(PetscSectionSetFieldName(section, one, 'v', ierr)) 77d8606c27SBarry Smith PetscCallA(PetscSectionSetFieldName(section, two, 'w', ierr)) 784820e4eaSBarry Smith if (size == 1) then 79d8606c27SBarry Smith PetscCallA(PetscSectionView(section, PETSC_VIEWER_STDOUT_WORLD, ierr)) 802e3d3ef9SMatthew G. Knepley end if 81c4762a1bSJed Brown! Tell the DM to use this data layout 82d8606c27SBarry Smith PetscCallA(DMSetLocalSection(dm, section, ierr)) 83c4762a1bSJed Brown! Create a Vec with this layout and view it 84d8606c27SBarry Smith PetscCallA(DMGetGlobalVector(dm, u, ierr)) 85d8606c27SBarry Smith PetscCallA(PetscViewerCreate(PETSC_COMM_WORLD, viewer, ierr)) 86d8606c27SBarry Smith PetscCallA(PetscViewerSetType(viewer, PETSCVIEWERVTK, ierr)) 87d8606c27SBarry Smith PetscCallA(PetscViewerFileSetName(viewer, 'sol.vtu', ierr)) 88d8606c27SBarry Smith PetscCallA(VecView(u, viewer, ierr)) 89d8606c27SBarry Smith PetscCallA(PetscViewerDestroy(viewer, ierr)) 90d8606c27SBarry Smith PetscCallA(DMRestoreGlobalVector(dm, u, ierr)) 91c4762a1bSJed Brown! Cleanup 92d8606c27SBarry Smith PetscCallA(PetscSectionDestroy(section, ierr)) 93d8606c27SBarry Smith PetscCallA(DMDestroy(dm, ierr)) 94c4762a1bSJed Brown 95d8606c27SBarry Smith PetscCallA(PetscFinalize(ierr)) 96c4762a1bSJed Brownend program DMPlexTestField 97c4762a1bSJed Brown 98c4762a1bSJed Brown!/*TEST 99*1c4a1a5bSMartin Diehl! 100c4762a1bSJed Brown! test: 101c4762a1bSJed Brown! suffix: 0 102c4762a1bSJed Brown! requires: triangle 10330602db0SMatthew G. Knepley! args: -info :~sys,mat: 104c4762a1bSJed Brown! 105c4762a1bSJed Brown! test: 1062e3d3ef9SMatthew G. Knepley! suffix: 0_2 1072e3d3ef9SMatthew G. Knepley! requires: triangle 1082e3d3ef9SMatthew G. Knepley! nsize: 2 1092e3d3ef9SMatthew G. Knepley! args: -info :~sys,mat,partitioner: 1102e3d3ef9SMatthew G. Knepley! 1112e3d3ef9SMatthew G. Knepley! test: 112c4762a1bSJed Brown! suffix: 1 113c4762a1bSJed Brown! requires: ctetgen 11430602db0SMatthew G. Knepley! args: -dm_plex_dim 3 -info :~sys,mat: 115c4762a1bSJed Brown! 116c4762a1bSJed Brown!TEST*/ 117