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