1 program main 2#include <petsc/finclude/petscdmplex.h> 3 use petscdmplex 4 implicit none 5! 6! 7 DM dm 8 PetscInt, dimension(4) :: EC 9 PetscInt, pointer :: pEC(:) 10 PetscInt, pointer :: pES(:) 11 PetscInt c, firstCell, numCells 12 PetscInt v, numVertices, numPoints 13 PetscInt i0,i4 14 PetscErrorCode ierr 15 16 i0 = 0 17 i4 = 4 18 19 PetscCallA(PetscInitialize(ierr)) 20 21 PetscCallA(DMPlexCreate(PETSC_COMM_WORLD, dm, ierr)) 22 firstCell = 0 23 numCells = 2 24 numVertices = 6 25 numPoints = numCells+numVertices 26 PetscCallA(DMPlexSetChart(dm, i0, numPoints, ierr)) 27 do c=firstCell,numCells-1 28 PetscCallA(DMPlexSetConeSize(dm, c, i4, ierr)) 29 end do 30 PetscCallA(DMSetUp(dm, ierr)) 31 32 EC(1) = 2 33 EC(2) = 3 34 EC(3) = 4 35 EC(4) = 5 36 c = 0 37 write(*,1000) 'cell EC 0',c,EC 38 1000 format (a,i4,50i4) 39 PetscCallA(DMPlexSetCone(dm, c , EC, ierr)) 40 PetscCallA(DMPlexGetCone(dm, c , pEC, ierr)) 41 write(*,1000) 'cell pEC 0',c,pEC 42 PetscCallA(DMPlexRestoreCone(dm, c, pEC, ierr)) 43 EC(1) = 4 44 EC(2) = 5 45 EC(3) = 6 46 EC(4) = 7 47 c = 1 48 write(*,1000) 'cell EC 1',c,EC 49 PetscCallA(DMPlexSetCone(dm, c , EC, ierr)) 50 PetscCallA(DMPlexGetCone(dm, c , pEC, ierr)) 51 write(*,1000) 'cell pEC 1',c,pEC 52 PetscCallA(DMPlexRestoreCone(dm, c, pEC, ierr)) 53 CHKMEMQ 54 55 PetscCallA(DMPlexSymmetrize(dm, ierr)) 56 PetscCallA(DMPlexStratify(dm, ierr)) 57 PetscCallA(DMPlexGetCone(dm, c , pEC, ierr)) 58 write(*,1000) 'cell pEC 3',c,pEC 59 PetscCallA(DMPlexRestoreCone(dm, c, pEC, ierr)) 60 61 v = 4 62 PetscCallA(DMPlexGetSupport(dm, v , pES, ierr)) 63 write(*,1000) 'vertex',v,pES 64 PetscCallA(DMPlexRestoreSupport(dm, v , pES, ierr)) 65 66 PetscCallA(DMDestroy(dm,ierr)) 67 PetscCallA(PetscFinalize(ierr)) 68 end 69 70! /*TEST 71! 72! test: 73! suffix: 0 74! 75! TEST*/ 76