Actual source code: ex1f90.F90

petsc-3.8.4 2018-03-24
Report Typos and Errors
  1:       program main
  2:  #include <petsc/finclude/petscdmplex.h>
  3:       use petscdmplex
  4:       implicit none
  5: !
  6: !
  7:       DM dm
  8:       PetscInt, target, 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

 16:       i0 = 0
 17:       i4 = 4

 19:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 20:      if (ierr .ne. 0) then
 21:         print*,'Unable to initialize PETSc'
 22:         stop
 23:       endif
 24:       call DMPlexCreate(PETSC_COMM_WORLD, dm, ierr);CHKERRA(ierr)
 25:       firstCell = 0
 26:       numCells = 2
 27:       numVertices = 6
 28:       numPoints = numCells+numVertices
 29:       call DMPlexSetChart(dm, i0, numPoints, ierr);CHKERRA(ierr)
 30:       do c=firstCell,numCells-1
 31:          call DMPlexSetConeSize(dm, c, i4, ierr);CHKERRA(ierr)
 32:       end do
 33:       call DMSetUp(dm, ierr);CHKERRA(ierr)

 35:       EC(1) = 2
 36:       EC(2) = 3
 37:       EC(3) = 4
 38:       EC(4) = 5
 39:       pEC => EC
 40:       c = 0
 41:       write(*,*) 'cell',c,pEC
 42:       call DMPlexSetCone(dm, c , pEC, ierr);CHKERRA(ierr)
 43:       call DMPlexGetCone(dm, c , pEC, ierr);CHKERRA(ierr)
 44:       write(*,*) 'cell',c,pEC
 45:       EC(1) = 4
 46:       EC(2) = 5
 47:       EC(3) = 6
 48:       EC(4) = 7
 49:       pEC => EC
 50:       c = 1
 51:       write(*,*) 'cell',c,pEC
 52:       call DMPlexSetCone(dm, c , pEC, ierr);CHKERRA(ierr)
 53:       call DMPlexGetCone(dm, c , pEC, ierr);CHKERRA(ierr)
 54:       write(*,*) 'cell',c,pEC
 55:       call DMPlexRestoreCone(dm, c , pEC, ierr);CHKERRA(ierr)

 57:       call DMPlexSymmetrize(dm, ierr);CHKERRA(ierr)
 58:       call DMPlexStratify(dm, ierr);CHKERRA(ierr)

 60:       v = 4
 61:       call DMPlexGetSupport(dm, v , pES, ierr);CHKERRA(ierr)
 62:       write(*,*) 'vertex',v,pES
 63:       call DMPlexRestoreSupport(dm, v , pES, ierr);CHKERRA(ierr)

 65:       call DMDestroy(dm,ierr);CHKERRA(ierr)
 66:       call PetscFinalize(ierr)
 67:       end

 69: ! /*TEST
 70: !
 71: ! test:
 72: !   suffix: 0
 73: !
 74: ! TEST*/