1#include <petsc/finclude/petscdmplex.h> 2program main 3 use petscdm 4 use petscdmplex 5 implicit none 6 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, nC 15 PetscInt i0, i1, i2, i3, i6, i7 16 PetscInt i8, i9, i10, i11 17 PetscErrorCode ierr 18 19 i0 = 0 20 i1 = 1 21 i2 = 2 22 i3 = 3 23 i6 = 6 24 i7 = 7 25 i8 = 8 26 i9 = 9 27 i10 = 10 28 i11 = 11 29 30 PetscCallA(PetscInitialize(ierr)) 31 32 PetscCallA(DMPlexCreate(PETSC_COMM_WORLD, dm, ierr)) 33 PetscCallA(PetscObjectSetName(dm, 'Mesh', ierr)) 34 dim = 2 35 PetscCallA(DMSetDimension(dm, dim, ierr)) 36 37! Make Doublet Mesh from Fig 2 of Flexible Representation of Computational Meshes, 38! except indexing is from 0 instead of 1 and we obey the new restrictions on 39! numbering: cells, vertices, faces, edges 40 PetscCallA(DMPlexSetChart(dm, i0, i11, ierr)) 41! cells 42 PetscCallA(DMPlexSetConeSize(dm, i0, i3, ierr)) 43 PetscCallA(DMPlexSetConeSize(dm, i1, i3, ierr)) 44! edges 45 PetscCallA(DMPlexSetConeSize(dm, i6, i2, ierr)) 46 PetscCallA(DMPlexSetConeSize(dm, i7, i2, ierr)) 47 PetscCallA(DMPlexSetConeSize(dm, i8, i2, ierr)) 48 PetscCallA(DMPlexSetConeSize(dm, i9, i2, ierr)) 49 PetscCallA(DMPlexSetConeSize(dm, i10, i2, ierr)) 50 51 PetscCallA(DMSetUp(dm, ierr)) 52 53 EC(1) = 6 54 EC(2) = 7 55 EC(3) = 8 56 pEC => EC 57 PetscCallA(DMPlexSetCone(dm, i0, pEC, ierr)) 58 EC(1) = 7 59 EC(2) = 9 60 EC(3) = 10 61 pEC => EC 62 PetscCallA(DMPlexSetCone(dm, i1, pEC, ierr)) 63 64 VE(1) = 2 65 VE(2) = 3 66 pVE => VE 67 PetscCallA(DMPlexSetCone(dm, i6, pVE, ierr)) 68 VE(1) = 3 69 VE(2) = 4 70 pVE => VE 71 PetscCallA(DMPlexSetCone(dm, i7, pVE, ierr)) 72 VE(1) = 4 73 VE(2) = 2 74 pVE => VE 75 PetscCallA(DMPlexSetCone(dm, i8, pVE, ierr)) 76 VE(1) = 3 77 VE(2) = 5 78 pVE => VE 79 PetscCallA(DMPlexSetCone(dm, i9, pVE, ierr)) 80 VE(1) = 5 81 VE(2) = 4 82 pVE => VE 83 PetscCallA(DMPlexSetCone(dm, i10, pVE, ierr)) 84 85 PetscCallA(DMPlexSymmetrize(dm, ierr)) 86 PetscCallA(DMPlexStratify(dm, ierr)) 87 PetscCallA(DMView(dm, PETSC_VIEWER_STDOUT_WORLD, ierr)) 88 89! Test Closure 90 do cell = 0, 1 91 PetscCallA(DMPlexGetTransitiveClosure(dm, cell, PETSC_TRUE, nC, nClosure, ierr)) 92! Different Fortran compilers print a different number of columns 93! per row producing different outputs in the test runs hence 94! do not print the nClosure 95 write (*, 1000) 'nClosure ', nClosure 961000 format(a, 30i4) 97 PetscCallA(DMPlexRestoreTransitiveClosure(dm, cell, PETSC_TRUE, nC, nClosure, ierr)) 98 end do 99 100! Test Join 101 size = 2 102 VE(1) = 6 103 VE(2) = 7 104 pVE => VE 105 PetscCallA(DMPlexGetJoin(dm, size, pVE, PETSC_NULL_INTEGER, nJoin, ierr)) 106 write (*, 1001) 'Join of', pVE 107 write (*, 1002) ' is', nJoin 108 PetscCallA(DMPlexRestoreJoin(dm, size, pVE, PETSC_NULL_INTEGER, nJoin, ierr)) 109 size = 2 110 VE(1) = 9 111 VE(2) = 7 112 pVE => VE 113 PetscCallA(DMPlexGetJoin(dm, size, pVE, PETSC_NULL_INTEGER, nJoin, ierr)) 114 write (*, 1001) 'Join of', pVE 1151001 format(a, 10i5) 116 write (*, 1002) ' is', nJoin 1171002 format(a, 10i5) 118 PetscCallA(DMPlexRestoreJoin(dm, size, pVE, PETSC_NULL_INTEGER, nJoin, ierr)) 119! Test Full Join 120 size = 3 121 EC(1) = 3 122 EC(2) = 4 123 EC(3) = 5 124 pEC => EC 125 PetscCallA(DMPlexGetFullJoin(dm, size, pEC, PETSC_NULL_INTEGER, nJoin, ierr)) 126 write (*, 1001) 'Full Join of', pEC 127 write (*, 1002) ' is', nJoin 128 PetscCallA(DMPlexRestoreJoin(dm, size, pEC, PETSC_NULL_INTEGER, nJoin, ierr)) 129! Test Meet 130 size = 2 131 VE(1) = 0 132 VE(2) = 1 133 pVE => VE 134 PetscCallA(DMPlexGetMeet(dm, size, pVE, PETSC_NULL_INTEGER, nMeet, ierr)) 135 write (*, 1001) 'Meet of', pVE 136 write (*, 1002) ' is', nMeet 137 PetscCallA(DMPlexRestoreMeet(dm, size, pVE, PETSC_NULL_INTEGER, nMeet, ierr)) 138 size = 2 139 VE(1) = 6 140 VE(2) = 7 141 pVE => VE 142 PetscCallA(DMPlexGetMeet(dm, size, pVE, PETSC_NULL_INTEGER, nMeet, ierr)) 143 write (*, 1001) 'Meet of', pVE 144 write (*, 1002) ' is', nMeet 145 PetscCallA(DMPlexRestoreMeet(dm, size, pVE, PETSC_NULL_INTEGER, nMeet, ierr)) 146 147 PetscCallA(DMDestroy(dm, ierr)) 148 PetscCallA(PetscFinalize(ierr)) 149end 150! 151!/*TEST 152! 153! test: 154! suffix: 0 155! 156!TEST*/ 157