<div dir="ltr"><div dir="ltr">On Thu, Mar 26, 2020 at 5:14 PM Yann Jobic <<a href="mailto:yann.jobic@univ-amu.fr" target="_blank">yann.jobic@univ-amu.fr</a>> wrote:<br></div><div class="gmail_quote"><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>
Le 3/23/2020 à 2:24 PM, Matthew Knepley a écrit :<br>
> On Wed, Mar 18, 2020 at 12:58 PM Yann Jobic <<a href="mailto:yann.jobic@univ-amu.fr" target="_blank">yann.jobic@univ-amu.fr</a> <br>
> <mailto:<a href="mailto:yann.jobic@univ-amu.fr" target="_blank">yann.jobic@univ-amu.fr</a>>> wrote:<br>
> <br>
>     Hi matt,<br>
> <br>
>     Le 3/17/2020 à 4:00 PM, Matthew Knepley a écrit :<br>
>      > On Mon, Mar 16, 2020 at 5:20 PM Yann Jobic<br>
>     <<a href="mailto:yann.jobic@univ-amu.fr" target="_blank">yann.jobic@univ-amu.fr</a> <mailto:<a href="mailto:yann.jobic@univ-amu.fr" target="_blank">yann.jobic@univ-amu.fr</a>><br>
>      > <mailto:<a href="mailto:yann.jobic@univ-amu.fr" target="_blank">yann.jobic@univ-amu.fr</a> <mailto:<a href="mailto:yann.jobic@univ-amu.fr" target="_blank">yann.jobic@univ-amu.fr</a>>>><br>
>     wrote:<br>
>      ><br>
>      >     Hi all,<br>
>      ><br>
>      >     I would like to implement a nodal DG with the DMPlex interface.<br>
>      >     Therefore, i must add the internal nodes to the DM (GLL<br>
>     nodes), with<br>
>      >     the<br>
>      >     constrains :<br>
>      >     1) Add them as solution points, with correct coordinates (and<br>
>     keep the<br>
>      >     good rotational ordering)<br>
>      >     2) Find the shared nodes at faces in order to compute the fluxes<br>
>      >     3) For parallel use, so synchronize the ghost node at each<br>
>     time steps<br>
>      ><br>
>      ><br>
>      > Let me get the fundamentals straight before advising, since I<br>
>     have never<br>
>      > implemented nodal DG.<br>
>      ><br>
>      >    1) What is shared?<br>
>     I need to duplicate an edge in 2D, or a facet in 3D, and to sync it<br>
>     after a time step, in order to compute the numerical fluxes<br>
>     (Lax-Friedrichs at the beginning).<br>
> <br>
> <br>
> I should have been more specific, but I think I see what you want. You <br>
> do not "share" unknowns between cells,<br>
> so all unknowns should be associated with some cell in the Section.<br>
> <br>
> You think of some cell unknowns as being "connected" to a face, so when <br>
> you want to calculate a flux, you need<br>
> the unknowns from the adjacent cell in order to do it. In order to do <br>
> this, I would partition with overlap=1, which<br>
> is what we do for finite volume, which has the same adjacency needs. You <br>
> might also set<br>
> <br>
> <a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DM/DMSetAdjacency.html" rel="noreferrer" target="_blank">https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DM/DMSetAdjacency.html</a><br>
> to PETSC_TRUE, PETSC_FALSE, but you are probably doing everything <br>
> matrix-free if you are using DG.<br>
> The above is optimal for FV, but not for DG because you communicate more <br>
> than you absolutely have to.<br>
> <br>
> A more complicated, but optimal, thing to do would be to assign interior <br>
> dofs to the cell, and two sets of dofs<br>
> to each face, one for each cell. Then you only communicate the face <br>
> dofs. Its just more bookkeeping for you,<br>
> but it will work in parallel just fine.<br>
I'm going this way.<br>
So i should use dm/impls/plex/examples/tutorials/ex1.c.html as reference.<br>
I should define the internal nodes on cells :<br>
I define the section with 3 fields (field 0 on cells, field 1 and 2 on <br>
faces), as :<br>
numComp[0] = Nr; /* Total number of dof per Cell */<br>
numDof[0*(dim+1)+0] = dim; /* defined over the Cell */<br>
And on the same section, the dofs at faces :<br>
numComp[1] = NumDofPerFace;<br>
numComp[2] = NumDofPerFace;<br>
numDof[1*(dim+1)+dim-1] = dim-1; /* internal dof of the cell */<br>
numDof[2*(dim+1)+dim-1] = dim-1; /* external dof of the cell */<br>
<br>
Is it a good way to create my section ?<br></blockquote><div><br></div><div>I would put them all in one field. A "field" is supposed to be a physical thing, like velocity or pressure.</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">
Thus, the data is duplicated for the faces, that means that i have to <br>
sync the internal Face dof at faces with their corresponding values from <br>
the internal one (at cells).<br></blockquote><div><br></div><div>That is not how I was envisioning it. Perhaps a drawing. Suppose you had DG2, then you have</div><div><br></div><div>        17              18</div><div>  7-----8-----9-----14----15</div><div>   |               |                |</div><div>16,3    4    5,6    11   12,13</div><div>   |               |                |</div><div>   1-----2-----3-----9-----10</div><div>          19            20</div><div><br></div><div>so each face gets 2 dofs, one for each cell.When doing a cell integral, you only use the dof that is for that cell.</div><div>The local-to-global would update the face dofs, so you would get each side.</div><div><br></div><div>There is a reordering when you extract the closure. I have written one for spectral elements. We would need</div><div>another here that ordered all the "other" face dofs to the end.</div><div><br></div><div>This seems a little complicated to me. Do you know how Andreas Klockner does it in Hedge? Or Tim Warburton?</div><div>I just want to make sure I am not missing an elegant way to handle this.</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">
Here field 1 is synchronised with field 0, locally.<br>
But the external Face dof, field 2, have to be synchronised with the <br>
values of the adjacent cell.<br>
Is it possible to use something like  DMPlexGetFaceFields ?<br>
Is there an example of such use of PetscSection and synchronisation <br>
process ?<br>
<br>
For the parallel part, should i use PetscSF object ?<br></blockquote><div><br></div><div>In parallel, integrals would be summed into the global vector, so each side has a 0 for the other face dof and the right contribution</div><div>for its face dof. Then both sides get both solution dofs. It seems to work in my head.</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 read your article "Mesh Algorithms for PDE with Sieve I: Mesh <br>
Distribution". But it's refereeing to Matthew G. Knepley and Dmitry A. <br>
Karpeev. Sieve implementation.<br>
Technical Report ANL/MCS to appear, Argonne National Laboratory,<br>
January 2008.<br>
I couldn't find it. It is freely available ?<br></blockquote><div><br></div><div>Don't bother reading that. There are later ones:</div><div><br></div><div> There are two pretty good sources:</div><div><br></div><div>  <a href="https://arxiv.org/abs/1505.04633" target="_blank">https://arxiv.org/abs/1505.04633</a></div><div>  <a href="https://arxiv.org/abs/1506.06194" target="_blank">https://arxiv.org/abs/1506.06194</a></div><div><br></div><div>The last one is a follow-on to this paper</div><div><br></div><div>  <a href="https://arxiv.org/abs/0908.4427" target="_blank">https://arxiv.org/abs/0908.4427</a></div><br class="gmail-Apple-interchange-newline"><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">
> <br>
> I don't think you need extra vertices, > or coordinates, and for output I<br>
> recommend using DMPlexProject() to get<br>
> the solution in some space that can be plotted like P1, or anything else <br>
> supported by your visualization.<br>
<br>
I would like to use DMplex as much as i can, as i would in the future <br>
refine locally the mesh.<br>
<br>
I hope you're good in this difficult situation (covid19),<br>
<br>
Best regards,<br>
<br>
Yann<br>
<br>
> <br>
>    Thanks,<br>
> <br>
>       Matt<br>
> <br>
>      ><br>
>      >        We have an implementation of spectral element ordering<br>
>      ><br>
>     (<a href="https://gitlab.com/petsc/petsc/-/blob/master/src/dm/impls/plex/examples/tutorials/ex6.c" rel="noreferrer" target="_blank">https://gitlab.com/petsc/petsc/-/blob/master/src/dm/impls/plex/examples/tutorials/ex6.c</a>).<br>
> <br>
>      > Those share<br>
>      >        the whole element boundary.<br>
>      ><br>
>      >    2) What ghosts do you need?<br>
>     In order to compute the numerical fluxes of one element, i need the<br>
>     values of the surrounding nodes connected to the adjacent elements.<br>
>      ><br>
>      >    3) You want to store real space coordinates for a quadrature?<br>
>     It should be basically the same as PetscFE of higher order.<br>
>     I add some vertex needed to compute a polynomal solution of the desired<br>
>     order. That means that if i have a N, order of the local approximation,<br>
>     i need 0.5*(N+1)*(N+2) vertex to store in the DMPlex (in 2D), in<br>
>     order to :<br>
>     1) have the correct number of dof<br>
>     2) use ghost nodes to sync the values of the vertex/edge/facet for<br>
>     1D/2D/3D problem<br>
>     2) save correctly the solution<br>
> <br>
>     Does it make sense to you ?<br>
> <br>
>     Maybe like<br>
>     <a href="https://www.mcs.anl.gov/petsc/petsc-current/src/ts/examples/tutorials/ex11.c.html" rel="noreferrer" target="_blank">https://www.mcs.anl.gov/petsc/petsc-current/src/ts/examples/tutorials/ex11.c.html</a><br>
>     With the use of the function SplitFaces, which i didn't fully<br>
>     understood<br>
>     so far.<br>
> <br>
>     Thanks,<br>
> <br>
>     Yann<br>
> <br>
>      ><br>
>      >        We usually define a quadrature on the reference element once.<br>
>      ><br>
>      >    Thanks,<br>
>      ><br>
>      >      Matt<br>
>      ><br>
>      >     I found elements of answers in those threads :<br>
>      ><br>
>     <a href="https://lists.mcs.anl.gov/pipermail/petsc-users/2016-August/029985.html" rel="noreferrer" target="_blank">https://lists.mcs.anl.gov/pipermail/petsc-users/2016-August/029985.html</a><br>
>      ><br>
>     <a href="https://lists.mcs.anl.gov/mailman/htdig/petsc-users/2019-October/039581.html" rel="noreferrer" target="_blank">https://lists.mcs.anl.gov/mailman/htdig/petsc-users/2019-October/039581.html</a><br>
>      ><br>
>      >     However, it's not clear for me where to begin.<br>
>      ><br>
>      >     Quoting Matt, i should :<br>
>      >     "  DMGetCoordinateDM(dm, &cdm);<br>
>      >         <Set field information into cdm instead of dm><br>
>      >        DMCreateLocalVector(cdm, &coordinatesLocal);<br>
>      >        <Fill in higher order coordinate values><br>
>      >        DMSetCoordinatesLocal(dm, coordinatesLocal);"<br>
>      ><br>
>      >     However, i will not create ghost nodes this way. And i'm not<br>
>     sure to<br>
>      >     keep the good ordering.<br>
>      >     This part should be implemented in the PetscFE interface, for<br>
>     high<br>
>      >     order<br>
>      >     discrete solutions.<br>
>      >     I did not succeed in finding the correct part of the source<br>
>     doing it.<br>
>      ><br>
>      >     Could you please give me some hint to begin correctly thoses<br>
>     tasks ?<br>
>      ><br>
>      >     Thanks,<br>
>      ><br>
>      >     Yann<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/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</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> <<a href="http://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a>><br>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr"><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>