[petsc-dev] Issues with DMComplexGetTransitiveClosure in Fortran
Chris Eldred
chris.eldred at gmail.com
Thu Aug 30 11:57:15 CDT 2012
Oh ok that makes sense. The nomenclature was confusing me.
On Thu, Aug 30, 2012 at 10:53 AM, Matthew Knepley <knepley at gmail.com> wrote:
> On Thu, Aug 30, 2012 at 11:49 AM, Chris Eldred <chris.eldred at gmail.com>
> wrote:
>>
>> That is what I have already been doing- so that is good! As far as the
>> ordering convention goes, for a 2D mesh I should have vertices, faces
>> (what I have been calling cells) and then edges? Ex. for the doublet
>> mesh it would be 0,1,2,3 (vertices); 4,5 (faces); 6,7,8,9,10 (edges)?
>
>
> No, I do it exactly as in ex2f90:
>
> cells: 0,1
> vertices: 2,3,4,5
> edges: 6,7,8,9,10
>
> Nomenclature:
>
> vertex: Complex in mesh of dimension 0
> cell: Complex in mesh of co-dimension 0 (so highest dimension
> possible)
> edge: Complex of dimension 1
> face: Complex of co-dimension 1 (notice it can also be an edge in 2D)
>
> Thanks,
>
> Matt
>
>>
>> On Wed, Aug 29, 2012 at 8:51 PM, Matthew Knepley <knepley at gmail.com>
>> wrote:
>> > On Wed, Aug 29, 2012 at 1:35 PM, Chris Eldred <chris.eldred at gmail.com>
>> > wrote:
>> >>
>> >> DMComplexGetTransitiveClosure (and Restore) work fine now. What are
>> >> the new restrictions on numbering mentioned in
>> >> src/dm/impls/complex/examples/tests/ex2f90.F?
>> >
>> >
>> > In the original Sieve, you could use any id for any point. This is very
>> > flexible, but
>> > hugely expensive in terms of memory, and somewhat expensive in terms of
>> > time.
>> >
>> > In DMComplex, each stratum in the DAG (it is technically a ranked poset)
>> > must
>> > be numbered contiguously. This means that vertices, faces, etc. are all
>> > numbered
>> > contiguously. It turns out that in practice, this is not an onerous
>> > requirement, and
>> > it allows very effective optimizations.
>> >
>> > In addition, we have the convention that we number cells, then vertices,
>> > then faces,
>> > then edges. Nothing in the interface demands it, but all of the PETSc
>> > code
>> > for FEM
>> > obeys this since it simplifies debugging.
>> >
>> > Matt
>> >
>> >>
>> >> On Tue, Aug 28, 2012 at 8:16 PM, Matthew Knepley <knepley at gmail.com>
>> >> wrote:
>> >> > On Tue, Aug 28, 2012 at 5:01 PM, Chris Eldred
>> >> > <chris.eldred at gmail.com>
>> >> > wrote:
>> >> >>
>> >> >> On Tue, Aug 28, 2012 at 3:49 PM, Matthew Knepley <knepley at gmail.com>
>> >> >> wrote:
>> >> >> > On Tue, Aug 28, 2012 at 4:39 PM, Chris Eldred
>> >> >> > <chris.eldred at gmail.com>
>> >> >> > wrote:
>> >> >> >>
>> >> >> >> DMComplexGetTransitiveClosure does not appear to have the calling
>> >> >> >> sequence (DM dm, PetscInt p, PetscBool useCone, PetscInt
>> >> >> >> *numPoints,
>> >> >> >> PetscInt *points[]) indicated in the documentation- instead it
>> >> >> >> has
>> >> >> >> DM
>> >> >> >> dm, PetscInt p, PetscBool useCone, PetscInt *points[]. In
>> >> >> >> addition,
>> >> >> >> it
>> >> >> >> does not appear to be returning the full closure- instead it just
>> >> >> >> gives cone + the point.
>> >> >> >
>> >> >> >
>> >> >> > 1) I will turn this code into a running example
>> >> >> >
>> >> >>
>> >> >> Awesome! I have been using the Doublet mesh to test DMComplex
>> >> >> routines
>> >> >> and my own functions, so that is great.
>> >> >
>> >> >
>> >> > Okay, there was a bug in the array handling. The test is
>> >> >
>> >> > src/dm/impls/complex/examples/tests/ex2f90.F
>> >> >
>> >> > Thanks,
>> >> >
>> >> > Matt
>> >> >
>> >> >>
>> >> >> > 2) I will fix the documentation. The interface change for Fortran
>> >> >> > comes
>> >> >> > because the size
>> >> >> > is embedded in the array pointer. If you think I should still
>> >> >> > return
>> >> >> > the
>> >> >> > size in Fortran,
>> >> >> > let me know.
>> >> >>
>> >> >> I can get the size from the array pointer so just returning that
>> >> >> works.
>> >> >>
>> >> >> >
>> >> >> > 3) I will also check to make sure the documentation is clear that
>> >> >> > this
>> >> >> > thing
>> >> >> > returns each
>> >> >> > point AND its orientation.
>> >> >> >
>> >> >> > Will write again when I have your thing running. It could be an
>> >> >> > error
>> >> >> > with
>> >> >> > the Fortran array
>> >> >> > processing I put in.
>> >> >> >
>> >> >> > Thanks,
>> >> >> >
>> >> >> > Matt
>> >> >> >
>> >> >> >>
>> >> >> >> Sample code is below:
>> >> >> >>
>> >> >> >> program main
>> >> >> >> implicit none
>> >> >> >> !
>> >> >> >> ! #echo '#include <finclude/petscsys.h>\n'
>> >> >> >> #echo '#include <finclude/petsc.h90>'
>> >> >> >> #echo '#include <finclude/petscdmcomplex.h>\n'
>> >> >> >> #echo '#include <finclude/petscdmcomplex.h90>'
>> >> >> >> ! #echo '#include <finclude/petscdm.h>\n'
>> >> >> >> ! #echo '#include <finclude/petscdm.h90>\n'
>> >> >> >>
>> >> >> >> !
>> >> >> >> DM dm
>> >> >> >> PetscInt, target, dimension(3) :: EC
>> >> >> >> PetscInt, target, dimension(2) :: VE
>> >> >> >> PetscInt, pointer :: pEC(:), pVE(:)
>> >> >> >> PetscInt, pointer :: nClosure(:)
>> >> >> >> PetscInt, pointer :: pClosure(:), pCone(:), pCone2(:)
>> >> >> >> PetscErrorCode ierr
>> >> >> >>
>> >> >> >> call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> call DMComplexCreate(PETSC_COMM_WORLD, dm, ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >>
>> >> >> >> !Make Doublet Mesh from Fig 2 of Flexible Representation of
>> >> >> >> Computational Meshes (except indexing is from 0 instead of 1)
>> >> >> >> call DMComplexSetChart(dm, 0, 11, ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >>
>> >> >> >> call DMComplexSetConeSize(dm, 9, 3,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> call DMComplexSetConeSize(dm, 10, 3,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> call DMComplexSetConeSize(dm, 4, 2,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> call DMComplexSetConeSize(dm, 5, 2,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> call DMComplexSetConeSize(dm, 6, 2,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> call DMComplexSetConeSize(dm, 7, 2,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> call DMComplexSetConeSize(dm, 8, 2,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >>
>> >> >> >> call DMSetUp(dm, ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >>
>> >> >> >> EC(1) = 4
>> >> >> >> EC(2) = 5
>> >> >> >> EC(3) = 6
>> >> >> >> pEC => EC
>> >> >> >> call DMComplexSetCone(dm, 9 , pEC, ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> EC(1) = 6
>> >> >> >> EC(2) = 7
>> >> >> >> EC(3) = 8
>> >> >> >> pEC => EC
>> >> >> >> call DMComplexSetCone(dm, 10 , pEC, ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >>
>> >> >> >> VE(1) = 0
>> >> >> >> VE(2) = 1
>> >> >> >> pVE => VE
>> >> >> >> call DMComplexSetCone(dm, 4 , pVE, ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> VE(1) = 0
>> >> >> >> VE(2) = 2
>> >> >> >> pVE => VE
>> >> >> >> call DMComplexSetCone(dm, 5 , pVE, ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> VE(1) = 1
>> >> >> >> VE(2) = 2
>> >> >> >> pVE => VE
>> >> >> >> call DMComplexSetCone(dm, 6 , pVE, ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> VE(1) = 2
>> >> >> >> VE(2) = 3
>> >> >> >> pVE => VE
>> >> >> >> call DMComplexSetCone(dm, 7 , pVE, ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> VE(1) = 1
>> >> >> >> VE(2) = 3
>> >> >> >> pVE => VE
>> >> >> >> call DMComplexSetCone(dm, 8 , pVE, ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >>
>> >> >> >> call DMComplexSymmetrize(dm,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> call DMComplexStratify(dm,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >>
>> >> >> >> !Test Stuff
>> >> >> >>
>> >> >> >> call
>> >> >> >> DMComplexGetTransitiveClosure(dm,10,PETSC_TRUE,nClosure,ierr)
>> >> >> >> !pClosure,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> write(*,*) nClosure
>> >> >> >> write(*,*) pClosure
>> >> >> >>
>> >> >> >> call DMDestroy(dm,ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> call PetscFinalize(ierr)
>> >> >> >> CHKERRQ(ierr)
>> >> >> >> end
>> >> >> >>
>> >> >> >>
>> >> >> >> The output is:
>> >> >> >> 10 0 6 0 7
>> >> >> >> 0 8
>> >> >> >>
>> >> >> >> and then there is an error in DMDestroy
>> >> >> >>
>> >> >> >> [0]PETSC ERROR: --------------------- Error Message
>> >> >> >> ------------------------------------
>> >> >> >> [0]PETSC ERROR: Object is in wrong state!
>> >> >> >> [0]PETSC ERROR: Work array still checked out!
>> >> >> >> [0]PETSC ERROR:
>> >> >> >>
>> >> >> >>
>> >> >> >>
>> >> >> >> ------------------------------------------------------------------------
>> >> >> >> [0]PETSC ERROR: Petsc Development HG revision:
>> >> >> >> b6fe82991deee4a0d3f9a20654bc7750a6b1fe0f HG Date: Mon Aug 27
>> >> >> >> 13:17:10
>> >> >> >> 2012 -0500
>> >> >> >> [0]PETSC ERROR: See docs/changes/index.html for recent updates.
>> >> >> >> [0]PETSC ERROR: See docs/faq.html for hints about trouble
>> >> >> >> shooting.
>> >> >> >> [0]PETSC ERROR: See docs/index.html for manual pages.
>> >> >> >> [0]PETSC ERROR:
>> >> >> >>
>> >> >> >>
>> >> >> >>
>> >> >> >> ------------------------------------------------------------------------
>> >> >> >> [0]PETSC ERROR: ./test on a arch-linu named Puget-101334 by user
>> >> >> >> Tue
>> >> >> >> Aug 28 15:31:06 2012
>> >> >> >> [0]PETSC ERROR: Libraries linked from
>> >> >> >> /home/user/Desktop/LIBRARIES/petsc-dev/arch-linux2-cxx-debug/lib
>> >> >> >> [0]PETSC ERROR: Configure run at Mon Aug 27 12:55:14 2012
>> >> >> >> [0]PETSC ERROR: Configure options --download-boost
>> >> >> >> --download-chaco
>> >> >> >> --download-ctetgen --download-f-blas-lapack --download-fiat
>> >> >> >> --download-generator --download-metis --download-ml
>> >> >> >> --download-mpich
>> >> >> >> --download-parmetis --download-scientificpython
>> >> >> >> --download-triangle
>> >> >> >> --with-clanguage=cxx --with-dynamic-loading
>> >> >> >> --with-shared-libraries
>> >> >> >> --with-sieve PETSC_ARCH=arch-linux2-cxx-debug
>> >> >> >> [0]PETSC ERROR:
>> >> >> >>
>> >> >> >>
>> >> >> >>
>> >> >> >> ------------------------------------------------------------------------
>> >> >> >> [0]PETSC ERROR: DMDestroy() line 234 in
>> >> >> >> /home/user/Desktop/LIBRARIES/petsc-dev/src/dm/interface/dm.c
>> >> >> >> application called MPI_Abort(MPI_COMM_WORLD, 73) - process 0
>> >> >> >> [unset]: aborting job:
>> >> >> >> application called MPI_Abort(MPI_COMM_WORLD, 73) - process 0
>> >> >> >>
>> >> >> >> Ideas?
>> >> >> >>
>> >> >> >> --
>> >> >> >> Chris Eldred
>> >> >> >> DOE Computational Science Graduate Fellow
>> >> >> >> Graduate Student, Atmospheric Science, Colorado State University
>> >> >> >> B.S. Applied Computational Physics, Carnegie Mellon University,
>> >> >> >> 2009
>> >> >> >> chris.eldred at gmail.com
>> >> >> >
>> >> >> >
>> >> >> >
>> >> >> >
>> >> >> > --
>> >> >> > What most experimenters take for granted before they begin their
>> >> >> > experiments
>> >> >> > is infinitely more interesting than any results to which their
>> >> >> > experiments
>> >> >> > lead.
>> >> >> > -- Norbert Wiener
>> >> >>
>> >> >>
>> >> >>
>> >> >> --
>> >> >> Chris Eldred
>> >> >> DOE Computational Science Graduate Fellow
>> >> >> Graduate Student, Atmospheric Science, Colorado State University
>> >> >> B.S. Applied Computational Physics, Carnegie Mellon University, 2009
>> >> >> chris.eldred at gmail.com
>> >> >
>> >> >
>> >> >
>> >> >
>> >> > --
>> >> > What most experimenters take for granted before they begin their
>> >> > experiments
>> >> > is infinitely more interesting than any results to which their
>> >> > experiments
>> >> > lead.
>> >> > -- Norbert Wiener
>> >>
>> >>
>> >>
>> >> --
>> >> Chris Eldred
>> >> DOE Computational Science Graduate Fellow
>> >> Graduate Student, Atmospheric Science, Colorado State University
>> >> B.S. Applied Computational Physics, Carnegie Mellon University, 2009
>> >> chris.eldred at gmail.com
>> >
>> >
>> >
>> >
>> > --
>> > What most experimenters take for granted before they begin their
>> > experiments
>> > is infinitely more interesting than any results to which their
>> > experiments
>> > lead.
>> > -- Norbert Wiener
>>
>>
>>
>> --
>> Chris Eldred
>> DOE Computational Science Graduate Fellow
>> Graduate Student, Atmospheric Science, Colorado State University
>> B.S. Applied Computational Physics, Carnegie Mellon University, 2009
>> chris.eldred at gmail.com
>
>
>
>
> --
> What most experimenters take for granted before they begin their experiments
> is infinitely more interesting than any results to which their experiments
> lead.
> -- Norbert Wiener
--
Chris Eldred
DOE Computational Science Graduate Fellow
Graduate Student, Atmospheric Science, Colorado State University
B.S. Applied Computational Physics, Carnegie Mellon University, 2009
chris.eldred at gmail.com
More information about the petsc-dev
mailing list