<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Wed, Mar 19, 2014 at 7:35 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">
  
    
  
  <div text="#000000" bgcolor="#FFFFFF">
    <br>
    <div>On 19/03/14 02:14, Matthew Knepley
      wrote:<br>
    </div>
    <blockquote type="cite">
      
      <div dir="ltr">
        <div class="gmail_extra">
          <div class="gmail_quote"><br>
            <blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">hi, I'm
              trying to use the DMPlexComputeCellGeometryFVM() function
              from a fortran program, but when I compile and link it
              complains:<br>
              <br>
              undefined reference to `dmplexcomputecellgeometryfvm_'<br>
              <br>
              I am running the 'next' branch of PETSc (having previously
              tried the master branch) and configured with
              --with-fortran-interfaces=1. At the top of my code I have
              #include <finclude/petsc.h90>.<br>
              <br>
              I noticed the fortran interfaces to some other functions
              (e.g. DMPlexCreateFromDAG) also seemed to be missing in
              older versions (e.g. the 'maint' branch)- is this one
              perhaps still to be added too?<br>
            </blockquote>
            <div><br>
            </div>
            <div>This is my fault. I will push it to 'next' today.</div>
          </div>
        </div>
      </div>
    </blockquote>
    <br>
    Thanks for doing that. I can access DMPlexComputeCellGeometryFVM()
    from fortran now.<br>
    <br>
    However at runtime it is giving me a segfault- perhaps I am not
    declaring the centroid and normal arrays properly? I also tried
    declaring them as ordinary fortran arrays instead of pointers, but
    that also segfaulted. It works OK when I call it from C. My fortran
    code is below- any clues much appreciated. The error it gives is:<br></div></blockquote><div><br></div><div>They need to be F90 pointers and they need to be associated with an allocated array. Take</div><div>a look at how I do it for the arguments to DMPlexCreateSection() in DMPlex ex1f90.F</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 text="#000000" bgcolor="#FFFFFF">

    [0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation
    Violation, probably memory access out of range<br>
    [0]PETSC ERROR: Try option -start_in_debugger or
    -on_error_attach_debugger<br>
    [0]PETSC ERROR: or see
    <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind</a>[0]PETSC
    ERROR: or try <a href="http://valgrind.org" target="_blank">http://valgrind.org</a> on G<br>
    NU/linux and Apple Mac OS X to find memory corruption
    errors                                                                    
    <br>
    [0]PETSC ERROR: likely location of problem given in stack below<br>
    [0]PETSC ERROR: ---------------------  Stack Frames
    ------------------------------------<br>
    [0]PETSC ERROR: Note: The EXACT line numbers in the stack are not
    available,<br>
    [0]PETSC ERROR:       INSTEAD the line number of the start of the
    function<br>
    [0]PETSC ERROR:       is given.<br>
    [0]PETSC ERROR: [0] F90Array1dAccess line 70
    /home/acro018/software/PETSc/code/src/sys/f90-src/f90_cwrap.c<br>
    [0]PETSC ERROR: --------------------- Error Message
    --------------------------------------------------------------<br>
    [0]PETSC ERROR: Signal received<br>
    [0]PETSC ERROR: See
    <a href="http://http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://http://www.mcs.anl.gov/petsc/documentation/faq.html</a> for
    trouble shooting.<br>
    [0]PETSC ERROR: Petsc Development GIT revision:
    v3.4.4-4938-g89a3bfc  GIT Date: 2014-03-19 17:00:11 +0000<br>
    [0]PETSC ERROR: testdag on a linux-gnu-c-opt named des108 by acro018
    Thu Mar 20 13:22:15 2014<br>
    [0]PETSC ERROR: Configure options --download-netcdf
    --download-exodusii --with-hdf5-dir=/usr/lib --with-numpy-dir=<br>
    [0]PETSC ERROR: #1 User provided function() line 0 in  unknown file<br>
    <br>
    - 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>
      PetscInt :: dim = 3<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>
      PetscReal ::  vol<br>
      PetscReal, pointer :: centroid(:), normal(:)<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>
    <br>
      call PetscInitialize(PETSC_NULL_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>
         call DMPlexComputeCellGeometryFVM(dm, i, vol, centroid, normal,
    ierr)<br>
         CHKERRQ(ierr)<br>
         write(*, '(a, i2, a, f8.4, a, 3(f8.4, 1x))'), &<br>
              'cell:', i, 'volume:', vol, 'centroid:', centroid(1),
    centroid(2), centroid(3)<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<span class="HOEnZb"><font color="#888888"><br>
    <br>
    <pre cols="72">-- 
Dr Adrian Croucher
Department of Engineering Science
University of Auckland
New Zealand
tel 64-9-373-7599 ext 84611
</pre>
  </font></span></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>