[petsc-users] DMPLEX project function error
yann JOBIC
yann.jobic at univ-amu.fr
Sun Oct 14 02:07:51 CDT 2018
Hi Ellen,
I've made an old code (i forgot in which version it's working) that
gives local coordinates and connectivity of a distributed DM.
I didn't check that it's working again. I'll also work on that on the
future.
If it can help, i send you my ugly code (i didn't clean it).
Yann
On 14/10/2018 06:36, Ellen M. Price wrote:
> Hi Matt,
>
> I ran into a problem implementing this. My integration code relies on
> knowing the points of each tetrahedron in the domain so I can turn an
> integral over all space into a sum of integrals over cells. I can't find
> an example of getting vertices for a given cell (I think TS ex11 comes
> closest, looks like it gets faces?). Is there an elegant way to do this,
> or an example I can follow?
>
> Ellen
>
>
> On 10/13/2018 04:26 AM, Matthew Knepley wrote:
>> On Fri, Oct 12, 2018 at 6:00 PM Ellen M. Price
>> <ellen.price at cfa.harvard.edu <mailto:ellen.price at cfa.harvard.edu>> wrote:
>>
>> Hi Matt,
>>
>> Thanks for the feedback! I'm always open to editorial comments. Do you
>> have a suggestion for a better way? I'm basically just trying to do an
>> integral over all space to get my boundary values for the FEM solve. I
>> don't know if PETSc has a way to do something like this?
>>
>>
>> So you are defining your BC by a boundary integral. I would just stick
>> that integration code
>> into a function and let that be called as your BC function. It will not
>> be called "extra" times,
>> just once for each dual vector, which in the case of P1 is every vertex.
>> This way it will work
>> for any element.
>>
>> Thanks,
>>
>> Matt
>>
>>
>> Ellen
>>
>>
>> On 10/12/2018 04:36 PM, Matthew Knepley wrote:
>> > On Fri, Oct 12, 2018 at 3:51 PM Ellen M. Price
>> > <ellen.price at cfa.harvard.edu <mailto:ellen.price at cfa.harvard.edu>
>> <mailto:ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>>> wrote:
>> >
>> > Hi Matt,
>> >
>> > Basically, I have a simple meshing script that generates a set
>> of points
>> > for my FEM mesh and outputs a C file that compiles into my
>> program.
>> > Along with those points is a static data field that I need to be
>> > associated with the correct data point. I don't have a
>> function that
>> > evaluates at each point, it's just pre-generated data.
>> >
>> >
>> > I cannot help but make editorial comments. Do not do this. It is
>> > incredibly fragile, since
>> > it depends on both a mesh and a discretization.
>> >
>> >
>> > Since I called DMPlexDistribute() early on, I think my mesh
>> points get
>> > re-ordered, right? So I need a way to re-order my data field
>> to match
>> > that ordering. Is there a good way to do this?
>> >
>> >
>> > However, if you want to do this, it is not hard. You do use
>> > DMPlexDistributeField
>> >
>> (https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DMPLEX/DMPlexDistributeField.html),
>> > but you need
>> > to be careful about the arguments. You want
>> >
>> > DMPlexDistribute(dmOriginal, 0, &migrationSF, &dmNew);
>> > DMPlexDistributeField(dmOriginal, migrationSF, originalSection,
>> > originalVec, newSection, newVec);
>> >
>> > Your Sections are easy to make, just put 1 dof on every vertex. The
>> > originalVec is just a Vec of your values
>> > (you can use VecCreateSeqWithArray() to wrap one around it), and the
>> > newVec you can get with
>> > DMCreateLocalVector(dmNew, &newVec) if you do DMSetSection(dmNew,
>> > newSection).
>> >
>> > Does this make sense?
>> >
>> > Thanks,
>> >
>> > Matt
>> >
>> >
>> > I tried using DMPlexDistributeField as a first attempt, but I
>> get some
>> > kind of error no matter what I do (presumably I'm not giving
>> the right
>> > input), and there aren't any examples in the docs.
>> >
>> > Ellen
>> >
>> >
>> > On 10/12/2018 11:31 AM, Matthew Knepley wrote:
>> > > On Fri, Oct 12, 2018 at 11:12 AM Ellen M. Price
>> > > <ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>
>> <mailto:ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>>
>> > <mailto:ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>
>> > <mailto:ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>>>> wrote:
>> > >
>> > > So I found *that* problem -- turns out that using
>> > DMPlexCreateSection
>> > > was messing things up, and it seems PETSc is handling
>> that part
>> > > internally?
>> > >
>> > >
>> > > DMPlexSection creates the Section based on the discretization
>> > > information in the DM (stored in a PetscDS).
>> > > If no discretization info is available, then it defaults to P0 I
>> > think,
>> > > which is what you got.
>> > >
>> > > After discretization info is put into the DM, if I ask for the
>> > Section,
>> > > the DM will create it automatically if it is missing.
>> > >
>> > >
>> > > However, I've run into a related problem now. I need to
>> > > project some field data onto the mesh; it's a static
>> array of
>> > > integration weights, and each one corresponds to a single
>> > point in the
>> > > mesh. Based on this answer:
>> > >
>> > >
>> >
>> https://lists.mcs.anl.gov/pipermail/petsc-users/2015-December/027964.html
>> > >
>> > > I think I need DMPlexDistributeField.
>> > >
>> > >
>> > > I don't think that is what you want. That is just sending data
>> > around.
>> > >
>> > > You want to define an FEM field. We use DMProject*() for
>> this. These
>> > > take as input either
>> > > a function of the coordinates (regular function) or another
>> FEM field.
>> > >
>> > > Can you tell me how your function is defined? I mean apart
>> from the
>> > > particular data you
>> > > have at points.
>> > >
>> > > Thanks,
>> > >
>> > > Matt
>> > >
>> > >
>> > > However, *that* function takes a
>> > > section as an argument, and I had to take out the part
>> where I
>> > create a
>> > > section to get things working! Calling
>> DMGetDefaultSection or
>> > > DMGetDefaultGlobalSection doesn't seem to help matters,
>> > either, because
>> > > then I seem to get a null section?
>> > >
>> > > So I either need to figure out how to set the section
>> properly
>> > or how to
>> > > get the section properly. As always, any help is
>> appreciated.
>> > >
>> > > Ellen
>> > >
>> > >
>> > > On 10/11/2018 05:41 PM, Matthew Knepley wrote:
>> > > > On Thu, Oct 11, 2018 at 4:39 PM Ellen M. Price
>> > > > <ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>
>> > <mailto:ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>>
>> > <mailto:ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>
>> > <mailto:ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>>>
>> > > <mailto:ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>
>> > <mailto:ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>>
>> > > <mailto:ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>
>> > <mailto:ellen.price at cfa.harvard.edu
>> <mailto:ellen.price at cfa.harvard.edu>>>>> wrote:
>> > > >
>> > > > I was working with a DMPLEX and FEM, following
>> SNES example
>> > > 12. I get
>> > > > the following error when I call DMProjectFunction,
>> but I
>> > don't
>> > > know what
>> > > > it means. Can anyone explain where I might have gone
>> > wrong, or
>> > > at least
>> > > > what this error is telling me? I think the point
>> closure
>> > size is
>> > > > correct, since my mesh is 3d simplex,
>> > > >
>> > > >
>> > > > Yes, if you have 3D SIMPLEX mesh and are using P1
>> elements,
>> > then you
>> > > > would have
>> > > > 4 dofs in the closure of a cell. The dual space
>> dimension is the
>> > > number
>> > > > of dual space
>> > > > basis vectors assigned to points in the closure. Since
>> it is
>> > 1, it
>> > > looks
>> > > > like you have a P0
>> > > > dual space. I assume you changed something in ex12?
>> > > >
>> > > > Thanks,
>> > > >
>> > > > Matt
>> > > >
>> > > >
>> > > > but what is the dual space
>> > > > dimension, and where might I have set it incorrectly?
>> > > >
>> > > > [0]PETSC ERROR: Nonconforming object sizes
>> > > > [0]PETSC ERROR: The section point closure size 4
>> != dual
>> > space
>> > > > dimension 1
>> > > > [0]PETSC ERROR: See
>> > > http://www.mcs.anl.gov/petsc/documentation/faq.html
>> > > > for trouble shooting.
>> > > > [0]PETSC ERROR: Petsc Release Version 3.9.2, May,
>> 20, 2018
>> > > > ...
>> > > > [0]PETSC ERROR: Configure options
>> > > > --prefix=/home/eprice/software/petsc-opt --with-hdf5=1
>> > > > --with-hdf5-dir=/home/eprice/software/hdf5-parallel
>> > --with-mpe=1
>> > > > --with-mpe-dir=/home/eprice/software/mpe
>> --with-debugging=0
>> > > > LDFLAGS="-pthread -lz" COPTFLAGS="-O3 -march=native
>> > -mtune=native"
>> > > > CXXOPTFLAGS="-O3 -march=native -mtune=native"
>> FOPTFLAGS="-O3
>> > > > -march=native -mtune=native" --with-mpi=1
>> > > > --with-mpi-dir=/home/eprice/software/mpich
>> --with-mumps=1
>> > > > --with-mumps-dir=/home/eprice/software/mumps
>> > --with-parmetis=1
>> > > > --with-parmetis-dir=/home/eprice/software/parmetis
>> > --with-metis=1
>> > > > --with-metis-dir=/home/eprice/software/parmetis
>> > --with-ptscotch=1
>> > > > --with-ptscotch-dir=/home/eprice/software/scotch
>> > > --with-scalapack=1
>> > > > --with-scalapack-dir=/home/eprice/software/scalapack
>> > > > [0]PETSC ERROR: #1 DMProjectLocal_Generic_Plex()
>> line 347 in
>> > > >
>> /h/sabriel0/src/petsc-3.9.2/src/dm/impls/plex/plexproject.c
>> > > > [0]PETSC ERROR: #2 DMProjectFunctionLocal_Plex()
>> line 428 in
>> > > >
>> /h/sabriel0/src/petsc-3.9.2/src/dm/impls/plex/plexproject.c
>> > > > [0]PETSC ERROR: #3 DMProjectFunctionLocal() line
>> 6265 in
>> > > > /h/sabriel0/src/petsc-3.9.2/src/dm/interface/dm.c
>> > > > [0]PETSC ERROR: #4 DMProjectFunction() line 6250 in
>> > > > /h/sabriel0/src/petsc-3.9.2/src/dm/interface/dm.c
>> > > > ...
>> > > >
>> > > > (I know this is an optimized PETSc build, but I
>> get the same
>> > > error from
>> > > > my debug build, it's just much slower.)
>> > > >
>> > > > Ellen
>> > > >
>> > > >
>> > > >
>> > > > --
>> > > > 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
>> > > >
>> > > > https://www.cse.buffalo.edu/~knepley/
>> <https://www.cse.buffalo.edu/%7Eknepley/>
>> > <https://www.cse.buffalo.edu/%7Eknepley/>
>> > > <https://www.cse.buffalo.edu/%7Eknepley/>
>> > > > <http://www.cse.buffalo.edu/%7Eknepley/>
>> > >
>> > >
>> > >
>> > > --
>> > > 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
>> > >
>> > > https://www.cse.buffalo.edu/~knepley/
>> <https://www.cse.buffalo.edu/%7Eknepley/>
>> > <https://www.cse.buffalo.edu/%7Eknepley/>
>> > > <http://www.cse.buffalo.edu/%7Eknepley/>
>> >
>> >
>> >
>> > --
>> > 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
>> >
>> > https://www.cse.buffalo.edu/~knepley/
>> <https://www.cse.buffalo.edu/%7Eknepley/>
>> > <http://www.cse.buffalo.edu/%7Eknepley/>
>>
>>
>>
>> --
>> 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
>>
>> https://www.cse.buffalo.edu/~knepley/
>> <http://www.cse.buffalo.edu/%7Eknepley/>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: dmload.c
Type: text/x-csrc
Size: 15415 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20181014/684c0179/attachment-0001.bin>
More information about the petsc-users
mailing list