[petsc-dev] indexing for periodic BC

[Barry Smith]

On Jun 19, 2011, at 9:46 AM, Matthew Knepley wrote:

> On Sun, Jun 19, 2011 at 1:38 PM, Jed Brown <jed at 59a2.org> wrote:
> On Fri, Jun 17, 2011 at 23:00, Barry Smith <bsmith at mcs.anl.gov> wrote:
>  It is possible that list of elements is not correct for periodic boundary conditions; I've never used it. I'm forwarding this to petsc-dev so Matt can check it and see if the elements make sense and what indices they routine (are the local, global? I don't know) if they are local then maybe you don't need to the MatSetValuesStencil() and just call MatSetValuesLocal().
> 
> The function returns local indices, so you would use MatSetValuesLocal(). I put a link to these in the documentation.
> 
> But the implementation is not correct for periodic domains. It should not be a great deal of work to implement. You would just need to teach that function about the "wrap-around" elements. But this is actually not trivial for a mapped grid. In no case can you just use the coordinates that were set on the DMDA, you have to extrapolate. But it's unclear how to extrapolate for a non-uniform mesh in order to make the periodic problem consistent. I don't think this is easy to rectify without changing the API for managing coordinates, but you can certainly make it work for regular Cartesian grids.
> 
> By the way, I think the whole premise of DMDAGetElements() is flawed. Shape functions are not provided so the user has to implement those themselves and manage consistency. The actual logic inside DMDAGetElements() is really small compared to defining quadratures and basis functions, so it's no particular effort for the user to manage it on their own. Additionally, managing it on their own will perform somewhat better (no intermediate indexing structures) and gives more flexibility in boundary conditions.
> 
> It might not be a bad idea to offer an easy-to-use, first order only, structured grid FEM API, but this ain't it.
> 
> Yes, I do not disagree that DMDAGetElements() is flawed. We just do not have anything better right now.

  For non-periodic domains it is trivial can just loop over the local mesh deciding on the elements. For periodic it is more complicated cause one has to decide where to process the ghost elements.

   This is a problem for Jungho and I right now because the VI problems we are solving should really have periodic boundary conditions. So we either fix the GetElements() (uniform grid is fine) or  she manages marching over the mesh herself.  What should it be?

   Barry


> 
>    Matt
> 
> -- 
> 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