Actual source code: ex2f90.F90

petsc-3.12.1 2019-10-22
Report Typos and Errors
  1:       program main
  2:  #include <petsc/finclude/petscdmplex.h>
  3:       use petscdmplex
  4:       use petscsys
  5:       implicit none

  7:       DM dm
  8:       PetscInt, target, dimension(3) :: EC
  9:       PetscInt, target, dimension(2) :: VE
 10:       PetscInt, pointer :: pEC(:), pVE(:)
 11:       PetscInt, pointer :: nClosure(:)
 12:       PetscInt, pointer :: nJoin(:)
 13:       PetscInt, pointer :: nMeet(:)
 14:       PetscInt       dim, cell, size
 15:       PetscInt i0,i1,i2,i3,i4,i5,i6,i7
 16:       PetscInt i8,i9,i10,i11
 17:       PetscErrorCode ierr

 19:       i0 = 0
 20:       i1 = 1
 21:       i2 = 2
 22:       i3 = 3
 23:       i4 = 4
 24:       i5 = 5
 25:       i6 = 6
 26:       i7 = 7
 27:       i8 = 8
 28:       i9 = 9
 29:       i10 = 10
 30:       i11 = 11

 32:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 33:       if (ierr .ne. 0) then
 34:         print*,'Unable to initialize PETSc'
 35:         stop
 36:       endif

 38:       call DMPlexCreate(PETSC_COMM_WORLD, dm, ierr);CHKERRA(ierr)
 39:       call PetscObjectSetName(dm, 'Mesh', ierr);CHKERRA(ierr)
 40:       dim = 2
 41:       call DMSetDimension(dm, dim, ierr);CHKERRA(ierr)

 43: ! Make Doublet Mesh from Fig 2 of Flexible Representation of Computational Meshes,
 44: ! except indexing is from 0 instead of 1 and we obey the new restrictions on
 45: ! numbering: cells, vertices, faces, edges
 46:       call DMPlexSetChart(dm, i0, i11, ierr);CHKERRA(ierr)
 47: !     cells
 48:       call DMPlexSetConeSize(dm, i0, i3, ierr);CHKERRA(ierr)
 49:       call DMPlexSetConeSize(dm, i1, i3, ierr);CHKERRA(ierr)
 50: !     edges
 51:       call DMPlexSetConeSize(dm,  i6, i2, ierr);CHKERRA(ierr)
 52:       call DMPlexSetConeSize(dm,  i7, i2, ierr);CHKERRA(ierr)
 53:       call DMPlexSetConeSize(dm,  i8, i2, ierr);CHKERRA(ierr)
 54:       call DMPlexSetConeSize(dm,  i9, i2, ierr);CHKERRA(ierr)
 55:       call DMPlexSetConeSize(dm, i10, i2, ierr);CHKERRA(ierr)

 57:       call DMSetUp(dm, ierr);CHKERRA(ierr)

 59:       EC(1) = 6
 60:       EC(2) = 7
 61:       EC(3) = 8
 62:       pEC => EC
 63:       call DMPlexSetCone(dm, i0, pEC, ierr);CHKERRA(ierr)
 64:       EC(1) = 7
 65:       EC(2) = 9
 66:       EC(3) = 10
 67:       pEC => EC
 68:       call DMPlexSetCone(dm, i1 , pEC, ierr);CHKERRA(ierr)

 70:       VE(1) = 2
 71:       VE(2) = 3
 72:       pVE => VE
 73:       call DMPlexSetCone(dm, i6 , pVE, ierr);CHKERRA(ierr)
 74:       VE(1) = 3
 75:       VE(2) = 4
 76:       pVE => VE
 77:       call DMPlexSetCone(dm, i7 , pVE, ierr);CHKERRA(ierr)
 78:       VE(1) = 4
 79:       VE(2) = 2
 80:       pVE => VE
 81:       call DMPlexSetCone(dm, i8 , pVE, ierr);CHKERRA(ierr)
 82:       VE(1) = 3
 83:       VE(2) = 5
 84:       pVE => VE
 85:       call DMPlexSetCone(dm, i9 , pVE, ierr);CHKERRA(ierr)
 86:       VE(1) = 5
 87:       VE(2) = 4
 88:       pVE => VE
 89:       call DMPlexSetCone(dm, i10 , pVE, ierr);CHKERRA(ierr)

 91:       call DMPlexSymmetrize(dm,ierr);CHKERRA(ierr)
 92:       call DMPlexStratify(dm,ierr);CHKERRA(ierr)
 93:       call DMView(dm, PETSC_VIEWER_STDOUT_WORLD, ierr);CHKERRA(ierr)

 95: !     Test Closure
 96:       do cell = 0,1
 97:          call DMPlexGetTransitiveClosure(dm,cell,PETSC_TRUE,nClosure,ierr);CHKERRA(ierr)
 98: !     Different Fortran compilers print a different number of columns
 99: !     per row producing different outputs in the test runs hence
100: !     do not print the nClosure
101: !       write(*,*) nClosure
102:        call DMPlexRestoreTransitiveClosure(dm,cell,PETSC_TRUE,nClosure,ierr);CHKERRA(ierr)
103:       end do

105: !     Test Join
106:       size  = 2
107:       VE(1) = 6
108:       VE(2) = 7
109:       pVE => VE
110:       call DMPlexGetJoin(dm, size, pVE, nJoin, ierr);CHKERRA(ierr)
111:       write(*,*) 'Join of',pVE,'is',nJoin
112:       call DMPlexRestoreJoin(dm, size, pVE, nJoin, ierr);CHKERRA(ierr)
113:       size  = 2
114:       VE(1) = 9
115:       VE(2) = 7
116:       pVE => VE
117:       call DMPlexGetJoin(dm, size, pVE, nJoin, ierr);CHKERRA(ierr)
118:       write(*,*) 'Join of',pVE,'is',nJoin
119:       call DMPlexRestoreJoin(dm, size, pVE, nJoin, ierr);CHKERRA(ierr)
120: !     Test Full Join
121:       size  = 3
122:       EC(1) = 3
123:       EC(2) = 4
124:       EC(3) = 5
125:       pEC => EC
126:       call DMPlexGetFullJoin(dm, size, pEC, nJoin, ierr);CHKERRA(ierr)
127:       write(*,*) 'Full Join of',pEC,'is',nJoin
128:       call DMPlexRestoreJoin(dm, size, pEC, nJoin, ierr);CHKERRA(ierr)
129: !     Test Meet
130:       size  = 2
131:       VE(1) = 0
132:       VE(2) = 1
133:       pVE => VE
134:       call DMPlexGetMeet(dm, size, pVE, nMeet, ierr);CHKERRA(ierr)
135:       write(*,*) 'Meet of',pVE,'is',nMeet
136:       call DMPlexRestoreMeet(dm, size, pVE, nMeet, ierr);CHKERRA(ierr)
137:       size  = 2
138:       VE(1) = 6
139:       VE(2) = 7
140:       pVE => VE
141:       call DMPlexGetMeet(dm, size, pVE, nMeet, ierr);CHKERRA(ierr)
142:       write(*,*) 'Meet of',pVE,'is',nMeet
143:       call DMPlexRestoreMeet(dm, size, pVE, nMeet, ierr);CHKERRA(ierr)

145:       call DMDestroy(dm, ierr);CHKERRA(ierr)
146:       call PetscFinalize(ierr)
147:       end
148: !
149: !/*TEST
150: !
151: !   test:
152: !     suffix: 0
153: !
154: !TEST*/