[petsc-users] Periodic boundary conditions in DMPlex
Pierre Seize
pierre.seize at onera.fr
Mon Oct 18 02:27:05 CDT 2021
On 15/10/21 22:30, Matthew Knepley wrote:
> On Fri, Oct 15, 2021 at 10:16 AM Pierre Seize <pierre.seize at onera.fr
> <mailto:pierre.seize at onera.fr>> wrote:
>
> I read everything again, I think I did not understand you at
> first. The first solution is to modify the DAG, so that the
> rightmost cell is linked to the leftmost face, right ? To do that,
> do I have to manually edit the DAG (the mesh is read from a file) ?
>
> Yes, the DAG would be modified if you want it for some particular mesh
> that we cannot read automatically. For example, we can read periodic
> GMsh meshes.
>
> If so, the mesh connectivity is like the one of a torus, then how
> does it work with the cells/faces coordinates ?
>
> You let the coordinate field be in a DG space, so that it can have jumps.
I'm not sure I fully understand this, what I'll do is that I will
experiment with a periodic GMSH mesh.
> Now I think the second method may be more straightforward. What's
> the idea ? Get the mapping with DMGetLocalToGlobalMapping, then
> create the mapping corresponding to the periodicity with
> ISLocalToGlobalMappingCreate, and finally
> ISLocalToGlobalMappingConcatenate ? I'm not sure this is the way,
> and I did not find something like DMSetLocalToGlobalMapping to
> restore the modified mapping.
>
> It is more complicated. We make the LocalToGlobalMap by looking at the
> PetscSection (essentially if gives function space information) and
> deciding which unknowns are removed from the global space.
> You would need to decide that unknowns constrained by periodicity are
> not present in the global space. Actually, this is not hard. You just
> mark them as constrained in the PetscSection, and all the layout
> functions will function correctly. However, then the LocalToGlobalMap
> will not be exactly right because the constrained unknowns will not be
> filled in (just like Dirichlet conditions). You would augment the
> map so that it fills those in by looking up their periodic
> counterparts. Jed has argued for this type of periodicity.
I also don't understand. One one side of my mesh I have : | ghost1 |
cell1 | cell2 | ... and on the other | cell_n-1 | cell_n | ghost_n |.
Are not the ghosts (from DMPlexConstructGhostCells) already constrained ?
I experimented on that too, I did:
DMGetSectionSF
PetscSFGetGraph
augment the graph by adding the local ghost cells to the leaves, and the
correct remote "true" cells to the roots
PetscSFSetGraph
and it seems to do what I want. Is this what you meant ? Is this a
correct way to use the PETSc objects ? Or is this just hacky and I'm
lucky it works ?
Pierre
>
> To me, the first kind is much more straightforward, but maybe this is
> because I find the topology code more clear.
>
> Thanks,
>
> Matt
>
> Pierre
>
>
> On 15/10/21 15:33, Pierre Seize wrote:
>>
>> When I first tried to handle the periodicity, I found the
>> DMPlexCreateBoxMesh function (I cannot find the cylinder one).
>>
>> From reading the sources, I understand that we do some work
>> either in DMPlexCreateCubeMesh_Internal or with DMSetPeriodicity.
>>
>> I tried to use DMSetPeriodicity before, for example with a 2x2
>> box on length 10. I did something like:
>>
>> const PetscReal maxCell[] = {2, 2};
>> const PetscReal L[] = {10, 10};
>> const DMBoundaryType bd[] = {DM_BOUNDARY_PERIODIC,
>> DM_BOUNDARY_PERIODIC};
>> DMSetPeriodicity(dm, PETSC_TRUE, maxCell, L, bd);
>> // or:
>> DMSetPeriodicity(dm, PETSC_TRUE, NULL, L, bd);
>>
>> but it did not work:
>>
>> VecSet(X, 1);
>> DMGetLocalVector(dm, &locX);
>> VecZeroEntries(locX);
>> DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX);
>> DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX);
>> VecView(locX, PETSC_VIEWER_STDOUT_WORLD);
>>
>> but the ghost cells values are all 0, only the real cells are 1.
>> So I guess DMSetPeriodicity alone is not sufficient to handle the
>> periodicity. Is there a way to do what I want ? That is set up my
>> DMPlex in a way that DMGlobalToLocalBegin/DMGlobalToLocalEnd do
>> exchange values between procs AND exchange the periodic values?
>>
>>
>> Thanks for the help
>>
>>
>> Pierre
>>
>>
>> On 15/10/21 14:03, Matthew Knepley wrote:
>>> On Fri, Oct 15, 2021 at 7:31 AM Pierre Seize
>>> <pierre.seize at onera.fr <mailto:pierre.seize at onera.fr>> wrote:
>>>
>>> It makes sense, thank you. In fact, both ways seems better
>>> than my way. The first one looks the most straightforward.
>>> Unfortunately I do not know how to implement either of them.
>>> Could you please direct me to the corresponding PETSc
>>> functions ?
>>>
>>> The first way is implemented for example in
>>> DMPlexCreateBoxMesh() and DMPlexCreateCylinderMesh(). The second
>>> is not implemented since
>>> there did not seem to be a general way to do it. I would help if
>>> you wanted to try coding it up.
>>>
>>> Thanks,
>>>
>>> Matt
>>>
>>> Pierre
>>>
>>>
>>> On 15/10/21 13:25, Matthew Knepley wrote:
>>>> On Fri, Oct 15, 2021 at 7:08 AM Pierre Seize
>>>> <pierre.seize at onera.fr <mailto:pierre.seize at onera.fr>> wrote:
>>>>
>>>> Hi,
>>>>
>>>> I'm writing a code using PETSc to solve NS equations
>>>> with FV on an
>>>> unstructured mesh. Therefore I use DMPlex.
>>>>
>>>> Regarding periodicity, I manage to implement it this way:
>>>>
>>>> - for each couple of boundaries that is linked with
>>>> periodicity, I
>>>> create a buffer vector with an ISLocalToGlobalMapping
>>>>
>>>> - then, when I need to fill the ghost cells
>>>> corresponding to the
>>>> periodicity, the i "true" cell of the local vector
>>>> fills the buffer
>>>> vector on location i with VecSetValuesBlockedLocal, then
>>>> VecAssemblyBegin/VecAssemblyEnd ensure each value is
>>>> send to the correct
>>>> location thanks to the mapping, then the i "ghost" cell
>>>> of the local
>>>> vector reads the vector on location i to get it's value.
>>>>
>>>>
>>>> It works, but it seems to me there is a better way,
>>>> with maybe PetscSF,
>>>> VecScatter, or something I don't know yet. Does anyone
>>>> have any advice ?
>>>>
>>>>
>>>> There are at least two other ways to handle this. First,
>>>> the method that is advocated in
>>>> Plex is to actually make a periodic geometry, meaning
>>>> connect the cells that are meant
>>>> to be connected. Then, if you partition with overlap = 1,
>>>> PetscGlobalToLocal() will fill in
>>>> these cell values automatically.
>>>>
>>>> Second, you could use a non-periodic geometry, but alter
>>>> the LocalToGlobal map such
>>>> that the cells gets filled in anyway. Many codes use this
>>>> scheme and it is straightforward
>>>> with Plex just by augmenting the map it makes automatically.
>>>>
>>>> Does this make sense?
>>>>
>>>> Thanks,
>>>>
>>>> Matt
>>>>
>>>> Pierre Seize
>>>>
>>>> --
>>>> 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/>
>>>
>>>
>>>
>>> --
>>> 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/>
>>
>
>
>
> --
> 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 --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20211018/53251752/attachment-0001.html>
More information about the petsc-users
mailing list