1*c4762a1bSJed Brown! setting up 3-D DMPlex using DMPlexCreateFromDAG(), DMPlexInterpolate() and 2*c4762a1bSJed Brown! DMPlexComputeCellGeometryFVM() 3*c4762a1bSJed Brown! Contributed by Adrian Croucher <a.croucher@auckland.ac.nz> 4*c4762a1bSJed Brown program main 5*c4762a1bSJed Brown use petscsys 6*c4762a1bSJed Brown use petscdmplex 7*c4762a1bSJed Brown#include <petsc/finclude/petscsys.h> 8*c4762a1bSJed Brown#include <petsc/finclude/petscdmplex.h> 9*c4762a1bSJed Brown implicit none 10*c4762a1bSJed Brown DM :: dm, dmi 11*c4762a1bSJed Brown PetscFV :: fvm 12*c4762a1bSJed Brown PetscInt, parameter :: dim = 3 13*c4762a1bSJed Brown 14*c4762a1bSJed Brown PetscInt :: depth = 1 15*c4762a1bSJed Brown PetscErrorCode :: ierr 16*c4762a1bSJed Brown PetscInt, dimension(2) :: numPoints 17*c4762a1bSJed Brown PetscInt, dimension(14) :: coneSize 18*c4762a1bSJed Brown PetscInt, dimension(16) :: cones 19*c4762a1bSJed Brown PetscInt, dimension(16) :: coneOrientations 20*c4762a1bSJed Brown PetscScalar, dimension(36) :: vertexCoords 21*c4762a1bSJed Brown PetscReal :: vol = 0. 22*c4762a1bSJed Brown PetscReal, target, dimension(dim) :: centroid 23*c4762a1bSJed Brown PetscReal, target, dimension(dim) :: normal 24*c4762a1bSJed Brown PetscReal, pointer :: pcentroid(:) 25*c4762a1bSJed Brown PetscReal, pointer :: pnormal(:) 26*c4762a1bSJed Brown 27*c4762a1bSJed Brown PetscReal, target, dimension(dim) :: v0 28*c4762a1bSJed Brown PetscReal, target, dimension(dim*dim) :: J 29*c4762a1bSJed Brown PetscReal, target, dimension(dim*dim) :: invJ 30*c4762a1bSJed Brown PetscReal, pointer :: pv0(:) 31*c4762a1bSJed Brown PetscReal, pointer :: pJ(:) 32*c4762a1bSJed Brown PetscReal, pointer :: pinvJ(:) 33*c4762a1bSJed Brown PetscReal :: detJ 34*c4762a1bSJed Brown 35*c4762a1bSJed Brown PetscInt :: i 36*c4762a1bSJed Brown 37*c4762a1bSJed Brown numPoints = [12, 2] 38*c4762a1bSJed Brown coneSize = [8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0] 39*c4762a1bSJed Brown cones = [2,5,4,3,6,7,8,9, 3,4,11,10,7,12,13,8] 40*c4762a1bSJed Brown coneOrientations = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,0] 41*c4762a1bSJed Brown vertexCoords = [ & 42*c4762a1bSJed Brown & -0.5,0.0,0.0, 0.0,0.0,0.0, 0.0,1.0,0.0, -0.5,1.0,0.0, & 43*c4762a1bSJed Brown & -0.5,0.0,1.0, 0.0,0.0,1.0, 0.0,1.0,1.0, -0.5,1.0,1.0, & 44*c4762a1bSJed Brown & 0.5,0.0,0.0, 0.5,1.0,0.0, 0.5,0.0,1.0, 0.5,1.0,1.0 ] 45*c4762a1bSJed Brown pcentroid => centroid 46*c4762a1bSJed Brown pnormal => normal 47*c4762a1bSJed Brown 48*c4762a1bSJed Brown pv0 => v0 49*c4762a1bSJed Brown pJ => J 50*c4762a1bSJed Brown pinvJ => invJ 51*c4762a1bSJed Brown 52*c4762a1bSJed Brown call PetscInitialize(PETSC_NULL_CHARACTER,ierr) 53*c4762a1bSJed Brown if (ierr .ne. 0) then 54*c4762a1bSJed Brown print*,'Unable to initialize PETSc' 55*c4762a1bSJed Brown stop 56*c4762a1bSJed Brown endif 57*c4762a1bSJed Brown 58*c4762a1bSJed Brown call DMPlexCreate(PETSC_COMM_WORLD, dm, ierr);CHKERRA(ierr) 59*c4762a1bSJed Brown call PetscObjectSetName(dm, 'testplex', ierr);CHKERRA(ierr) 60*c4762a1bSJed Brown call DMSetDimension(dm, dim, ierr);CHKERRA(ierr) 61*c4762a1bSJed Brown 62*c4762a1bSJed Brown call DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones,coneOrientations, vertexCoords, ierr);CHKERRA(ierr) 63*c4762a1bSJed Brown 64*c4762a1bSJed Brown call DMPlexInterpolate(dm, dmi, ierr);CHKERRA(ierr) 65*c4762a1bSJed Brown call DMPlexCopyCoordinates(dm, dmi, ierr);CHKERRA(ierr) 66*c4762a1bSJed Brown call DMDestroy(dm, ierr);CHKERRA(ierr) 67*c4762a1bSJed Brown dm = dmi 68*c4762a1bSJed Brown 69*c4762a1bSJed Brown call DMView(dm, PETSC_VIEWER_STDOUT_WORLD, ierr);CHKERRA(ierr) 70*c4762a1bSJed Brown 71*c4762a1bSJed Brown do i = 0, 1 72*c4762a1bSJed Brown call DMPlexComputeCellGeometryFVM(dm, i, vol, pcentroid, pnormal, ierr);CHKERRA(ierr) 73*c4762a1bSJed Brown write(*, '(a, i2, a, f8.4, a, 3(f8.4, 1x))') 'cell: ', i, ' volume: ', vol, ' centroid: ',pcentroid(1), pcentroid(2), pcentroid(3) 74*c4762a1bSJed Brown call DMPlexComputeCellGeometryAffineFEM(dm, i, pv0, pJ, pinvJ,detJ, ierr);CHKERRA(ierr) 75*c4762a1bSJed Brown end do 76*c4762a1bSJed Brown 77*c4762a1bSJed Brown call PetscFVCreate(PETSC_COMM_WORLD, fvm, ierr);CHKERRA(ierr) 78*c4762a1bSJed Brown call PetscFVSetUp(fvm, ierr);CHKERRA(ierr) 79*c4762a1bSJed Brown call PetscFVDestroy(fvm, ierr);CHKERRA(ierr) 80*c4762a1bSJed Brown 81*c4762a1bSJed Brown call DMDestroy(dm, ierr);CHKERRA(ierr) 82*c4762a1bSJed Brown call PetscFinalize(ierr) 83*c4762a1bSJed Brown end program main 84*c4762a1bSJed Brown 85*c4762a1bSJed Brown!/*TEST 86*c4762a1bSJed Brown! 87*c4762a1bSJed Brown! test: 88*c4762a1bSJed Brown! suffix: 0 89*c4762a1bSJed Brown! 90*c4762a1bSJed Brown!TEST*/ 91