<div dir="ltr"><div dir="ltr"><div dir="ltr"><div class="gmail_quote"><div dir="ltr">On Fri, Oct 12, 2018 at 3:51 PM Ellen M. Price <<a href="mailto:ellen.price@cfa.harvard.edu">ellen.price@cfa.harvard.edu</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">Hi Matt,<br>
<br>
Basically, I have a simple meshing script that generates a set of points<br>
for my FEM mesh and outputs a C file that compiles into my program.<br>
Along with those points is a static data field that I need to be<br>
associated with the correct data point. I don't have a function that<br>
evaluates at each point, it's just pre-generated data.<br></blockquote><div><br></div><div>I cannot help but make editorial comments. Do not do this. It is incredibly fragile, since</div><div>it depends on both a mesh and a discretization.</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
Since I called DMPlexDistribute() early on, I think my mesh points get<br>
re-ordered, right? So I need a way to re-order my data field to match<br>
that ordering. Is there a good way to do this?<br></blockquote><div><br></div><div>However, if you want to do this, it is not hard. You do use DMPlexDistributeField (<a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DMPLEX/DMPlexDistributeField.html">https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DMPLEX/DMPlexDistributeField.html</a>), but you need</div><div>to be careful about the arguments. You want</div><div><br></div><div> DMPlexDistribute(dmOriginal, 0, &migrationSF, &dmNew);</div><div> DMPlexDistributeField(dmOriginal, migrationSF, originalSection, originalVec, newSection, newVec);</div><div><br></div><div>Your Sections are easy to make, just put 1 dof on every vertex. The originalVec is just a Vec of your values</div><div>(you can use VecCreateSeqWithArray() to wrap one around it), and the newVec you can get with</div><div>DMCreateLocalVector(dmNew, &newVec) if you do DMSetSection(dmNew, newSection).</div><div><br></div><div>Does this make sense?</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
I tried using DMPlexDistributeField as a first attempt, but I get some<br>
kind of error no matter what I do (presumably I'm not giving the right<br>
input), and there aren't any examples in the docs.<br>
<br>
Ellen<br>
<br>
<br>
On 10/12/2018 11:31 AM, Matthew Knepley wrote:<br>
> On Fri, Oct 12, 2018 at 11:12 AM Ellen M. Price<br>
> <<a href="mailto:ellen.price@cfa.harvard.edu" target="_blank">ellen.price@cfa.harvard.edu</a> <mailto:<a href="mailto:ellen.price@cfa.harvard.edu" target="_blank">ellen.price@cfa.harvard.edu</a>>> wrote:<br>
> <br>
> So I found *that* problem -- turns out that using DMPlexCreateSection<br>
> was messing things up, and it seems PETSc is handling that part<br>
> internally? <br>
> <br>
> <br>
> DMPlexSection creates the Section based on the discretization<br>
> information in the DM (stored in a PetscDS).<br>
> If no discretization info is available, then it defaults to P0 I think,<br>
> which is what you got.<br>
> <br>
> After discretization info is put into the DM, if I ask for the Section,<br>
> the DM will create it automatically if it is missing.<br>
> <br>
> <br>
> However, I've run into a related problem now. I need to<br>
> project some field data onto the mesh; it's a static array of<br>
> integration weights, and each one corresponds to a single point in the<br>
> mesh. Based on this answer:<br>
> <br>
> <a href="https://lists.mcs.anl.gov/pipermail/petsc-users/2015-December/027964.html" rel="noreferrer" target="_blank">https://lists.mcs.anl.gov/pipermail/petsc-users/2015-December/027964.html</a><br>
> <br>
> I think I need DMPlexDistributeField.<br>
> <br>
> <br>
> I don't think that is what you want. That is just sending data around.<br>
> <br>
> You want to define an FEM field. We use DMProject*() for this. These<br>
> take as input either<br>
> a function of the coordinates (regular function) or another FEM field.<br>
> <br>
> Can you tell me how your function is defined? I mean apart from the<br>
> particular data you<br>
> have at points.<br>
> <br>
> Thanks,<br>
> <br>
> Matt<br>
> <br>
> <br>
> However, *that* function takes a<br>
> section as an argument, and I had to take out the part where I create a<br>
> section to get things working! Calling DMGetDefaultSection or<br>
> DMGetDefaultGlobalSection doesn't seem to help matters, either, because<br>
> then I seem to get a null section?<br>
> <br>
> So I either need to figure out how to set the section properly or how to<br>
> get the section properly. As always, any help is appreciated.<br>
> <br>
> Ellen<br>
> <br>
> <br>
> On 10/11/2018 05:41 PM, Matthew Knepley wrote:<br>
> > On Thu, Oct 11, 2018 at 4:39 PM Ellen M. Price<br>
> > <<a href="mailto:ellen.price@cfa.harvard.edu" target="_blank">ellen.price@cfa.harvard.edu</a> <mailto:<a href="mailto:ellen.price@cfa.harvard.edu" target="_blank">ellen.price@cfa.harvard.edu</a>><br>
> <mailto:<a href="mailto:ellen.price@cfa.harvard.edu" target="_blank">ellen.price@cfa.harvard.edu</a><br>
> <mailto:<a href="mailto:ellen.price@cfa.harvard.edu" target="_blank">ellen.price@cfa.harvard.edu</a>>>> wrote:<br>
> ><br>
> > I was working with a DMPLEX and FEM, following SNES example<br>
> 12. I get<br>
> > the following error when I call DMProjectFunction, but I don't<br>
> know what<br>
> > it means. Can anyone explain where I might have gone wrong, or<br>
> at least<br>
> > what this error is telling me? I think the point closure size is<br>
> > correct, since my mesh is 3d simplex,<br>
> ><br>
> ><br>
> > Yes, if you have 3D SIMPLEX mesh and are using P1 elements, then you<br>
> > would have<br>
> > 4 dofs in the closure of a cell. The dual space dimension is the<br>
> number<br>
> > of dual space<br>
> > basis vectors assigned to points in the closure. Since it is 1, it<br>
> looks<br>
> > like you have a P0<br>
> > dual space. I assume you changed something in ex12?<br>
> ><br>
> > Thanks,<br>
> ><br>
> > Matt<br>
> > <br>
> ><br>
> > but what is the dual space<br>
> > dimension, and where might I have set it incorrectly?<br>
> ><br>
> > [0]PETSC ERROR: Nonconforming object sizes<br>
> > [0]PETSC ERROR: The section point closure size 4 != dual space<br>
> > dimension 1<br>
> > [0]PETSC ERROR: See<br>
> <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" rel="noreferrer" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html</a><br>
> > for trouble shooting.<br>
> > [0]PETSC ERROR: Petsc Release Version 3.9.2, May, 20, 2018<br>
> > ...<br>
> > [0]PETSC ERROR: Configure options<br>
> > --prefix=/home/eprice/software/petsc-opt --with-hdf5=1<br>
> > --with-hdf5-dir=/home/eprice/software/hdf5-parallel --with-mpe=1<br>
> > --with-mpe-dir=/home/eprice/software/mpe --with-debugging=0<br>
> > LDFLAGS="-pthread -lz" COPTFLAGS="-O3 -march=native -mtune=native"<br>
> > CXXOPTFLAGS="-O3 -march=native -mtune=native" FOPTFLAGS="-O3<br>
> > -march=native -mtune=native" --with-mpi=1<br>
> > --with-mpi-dir=/home/eprice/software/mpich --with-mumps=1<br>
> > --with-mumps-dir=/home/eprice/software/mumps --with-parmetis=1<br>
> > --with-parmetis-dir=/home/eprice/software/parmetis --with-metis=1<br>
> > --with-metis-dir=/home/eprice/software/parmetis --with-ptscotch=1<br>
> > --with-ptscotch-dir=/home/eprice/software/scotch<br>
> --with-scalapack=1<br>
> > --with-scalapack-dir=/home/eprice/software/scalapack<br>
> > [0]PETSC ERROR: #1 DMProjectLocal_Generic_Plex() line 347 in<br>
> > /h/sabriel0/src/petsc-3.9.2/src/dm/impls/plex/plexproject.c<br>
> > [0]PETSC ERROR: #2 DMProjectFunctionLocal_Plex() line 428 in<br>
> > /h/sabriel0/src/petsc-3.9.2/src/dm/impls/plex/plexproject.c<br>
> > [0]PETSC ERROR: #3 DMProjectFunctionLocal() line 6265 in<br>
> > /h/sabriel0/src/petsc-3.9.2/src/dm/interface/dm.c<br>
> > [0]PETSC ERROR: #4 DMProjectFunction() line 6250 in<br>
> > /h/sabriel0/src/petsc-3.9.2/src/dm/interface/dm.c<br>
> > ...<br>
> ><br>
> > (I know this is an optimized PETSc build, but I get the same<br>
> error from<br>
> > my debug build, it's just much slower.)<br>
> ><br>
> > Ellen<br>
> ><br>
> ><br>
> ><br>
> > --<br>
> > What most experimenters take for granted before they begin their<br>
> > experiments is infinitely more interesting than any results to which<br>
> > their experiments lead.<br>
> > -- Norbert Wiener<br>
> ><br>
> > <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br>
> <<a href="https://www.cse.buffalo.edu/%7Eknepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/%7Eknepley/</a>><br>
> > <<a href="http://www.cse.buffalo.edu/%7Eknepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/%7Eknepley/</a>><br>
> <br>
> <br>
> <br>
> -- <br>
> What most experimenters take for granted before they begin their<br>
> experiments is infinitely more interesting than any results to which<br>
> their experiments lead.<br>
> -- Norbert Wiener<br>
> <br>
> <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br>
> <<a href="http://www.cse.buffalo.edu/%7Eknepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/%7Eknepley/</a>><br>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>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><br></div><div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div></div></div>