<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Thu, Mar 20, 2014 at 9:04 PM, Adrian Croucher <span dir="ltr"><<a href="mailto:a.croucher@auckland.ac.nz" target="_blank">a.croucher@auckland.ac.nz</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">hi,<div class=""><br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
They need to be F90 pointers and they need to be associated with an<br>
allocated array. Take<br>
a look at how I do it for the arguments to DMPlexCreateSection() in DMPlex<br>
ex1f90.F<br>
</blockquote></div>
Thanks, I tried that and it nearly works... it gives the correct volumes, but it gives centroids all zero, which isn't right. My amended code is below.<br>
<br>
If I debug into DMPlexComputeGeometryFVM_3D_<u></u>Internal() using gdb, and try 'p *centroid', it says 'Cannot access memory at address 0x0', which looks unhealthy.<br>
<br>
There is something a bit weird going on, because if I don't initialize the centroid and normal arrays, or even if I initialize them in the main code rather than in the variable declaration, I get a segfault in DMPlexComputeGeometryFVM_3D_<u></u>Internal().<br>

<br>
Guess I'm still doing something subtly wrong here?</blockquote><div><br></div><div>The interface needed to be declared. I pushed the fix, and your example to PETSc (ex3f90).</div><div><br></div><div>  Thanks,</div><div>
<br></div><div>      Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div class="">
- Adrian<br>
<br>
--<br>
program main<br>
<br>
! setting up 3-D DMPlex using DMPlexCreateFromDAG(), DMPlexInterpolate() and<br>
! DMPlexComputeCellGeometryFVM()<br>
<br>
  implicit none<br>
<br>
#include <finclude/petsc.h90><br>
<br>
  DM :: dm, dmi<br></div>
  PetscInt, parameter :: dim = 3<div class=""><br>
  PetscInt :: depth = 1<br>
  PetscErrorCode :: ierr<br>
  PetscInt, allocatable :: numPoints(:)<br>
  PetscInt, allocatable :: coneSize(:)<br>
  PetscInt, allocatable :: cones(:)<br>
  PetscInt, allocatable :: coneOrientations(:)<br>
  PetscScalar, allocatable :: vertexCoords(:)<br></div>
  PetscReal ::  vol = 0.<br>
  PetscReal, target, dimension(dim)  :: centroid = 0., normal = 0.<br>
  PetscReal, pointer :: pcentroid(:), pnormal(:)<div class=""><br>
  PetscInt :: i<br>
<br>
  numPoints = [12, 2]<br>
  coneSize  = [8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]<br>
  cones = [2,5,4,3,6,7,8,9,  3,4,11,10,7,12,13,8]<br>
  coneOrientations = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]<br>
  vertexCoords = [ &<br>
       -0.5,0.0,0.0, 0.0,0.0,0.0, 0.0,1.0,0.0, -0.5,1.0,0.0, &<br>
       -0.5,0.0,1.0, 0.0,0.0,1.0, 0.0,1.0,1.0, -0.5,1.0,1.0, &<br>
        0.5,0.0,0.0, 0.5,1.0,0.0, 0.5,0.0,1.0,  0.5,1.0,1.0 ]<br></div>
  pcentroid => centroid<br>
  pnormal => normal<div class=""><br>
<br>
  call PetscInitialize(PETSC_NULL_<u></u>CHARACTER,ierr); CHKERRQ(ierr)<br>
<br>
  call DMPlexCreate(PETSC_COMM_WORLD, dm, ierr); CHKERRQ(ierr)<br>
  call PetscObjectSetName(dm, "testplex", ierr); CHKERRQ(ierr)<br>
  call DMPlexSetDimension(dm, dim, ierr); CHKERRQ(ierr)<br>
<br>
  call DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, &<br>
       coneOrientations, vertexCoords, ierr); CHKERRQ(ierr)<br>
<br>
  call DMPlexInterpolate(dm, dmi, ierr); CHKERRQ(ierr)<br>
  call DMPlexCopyCoordinates(dm, dmi, ierr); CHKERRQ(ierr)<br>
  call DMDestroy(dm, ierr); CHKERRQ(ierr)<br>
  dm = dmi<br>
<br>
  call DMView(dm, PETSC_VIEWER_STDOUT_WORLD, ierr); CHKERRQ(ierr)<br>
<br>
  do i = 0, 1<br></div>
     call DMPlexComputeCellGeometryFVM(<u></u>dm, i, vol, pcentroid, pnormal, ierr)<div class=""><br>
     CHKERRQ(ierr)<br>
     write(*, '(a, i2, a, f8.4, a, 3(f8.4, 1x))'), &<br></div>
          'cell: ', i, ' volume: ', vol, ' centroid: ', &<br>
          pcentroid(1), pcentroid(2), pcentroid(3)<div class="HOEnZb"><div class="h5"><br>
  end do<br>
<br>
  call DMDestroy(dm, ierr); CHKERRQ(ierr)<br>
  call PetscFinalize(ierr); CHKERRQ(ierr)<br>
  deallocate(numPoints, coneSize, cones, coneOrientations, vertexCoords)<br>
<br>
end program main<br>
<br>
-- <br>
Dr Adrian Croucher<br>
Department of Engineering Science<br>
University of Auckland<br>
New Zealand<br>
tel 64-9-373-7599 ext 84611<br>
<br>
</div></div></blockquote></div><br><br clear="all"><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener
</div></div>