Actual source code: ex1f90.F
petsc-3.7.5 2017-01-01
1: program main
2: implicit none
3: !
4: #include <petsc/finclude/petsc.h90>
5: !
6: DM dm
7: PetscInt, target, dimension(4) :: EC
8: PetscInt, pointer :: pEC(:)
9: PetscInt, pointer :: pES(:)
10: PetscInt c, firstCell, numCells
11: PetscInt v, numVertices, numPoints
12: PetscInt i0,i4
13: PetscErrorCode ierr
15: i0 = 0
16: i4 = 4
18: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
19: call DMPlexCreate(PETSC_COMM_WORLD, dm, ierr)
20: firstCell = 0
21: numCells = 2
22: numVertices = 6
23: numPoints = numCells+numVertices
24: call DMPlexSetChart(dm, i0, numPoints, ierr)
25: do c=firstCell,numCells-1
26: call DMPlexSetConeSize(dm, c, i4, ierr)
27: end do
28: call DMSetUp(dm, ierr)
30: EC(1) = 2
31: EC(2) = 3
32: EC(3) = 4
33: EC(4) = 5
34: pEC => EC
35: c = 0
36: write(*,*) 'cell',c,pEC
37: call DMPlexSetCone(dm, c , pEC, ierr)
38: CHKERRQ(ierr)
39: call DMPlexGetCone(dm, c , pEC, ierr)
40: CHKERRQ(ierr)
41: write(*,*) 'cell',c,pEC
42: EC(1) = 4
43: EC(2) = 5
44: EC(3) = 6
45: EC(4) = 7
46: pEC => EC
47: c = 1
48: write(*,*) 'cell',c,pEC
49: call DMPlexSetCone(dm, c , pEC, ierr)
50: CHKERRQ(ierr)
51: call DMPlexGetCone(dm, c , pEC, ierr)
52: CHKERRQ(ierr)
53: write(*,*) 'cell',c,pEC
54: call DMPlexRestoreCone(dm, c , pEC, ierr)
55: CHKERRQ(ierr)
57: call DMPlexSymmetrize(dm, ierr)
58: call DMPlexStratify(dm, ierr)
60: v = 4
61: call DMPlexGetSupport(dm, v , pES, ierr)
62: CHKERRQ(ierr)
63: write(*,*) 'vertex',v,pES
64: call DMPlexRestoreSupport(dm, v , pES, ierr)
65: CHKERRQ(ierr)
67: call DMDestroy(dm,ierr)
68: call PetscFinalize(ierr)
69: end