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