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