On 15/10/21 22:30, Matthew Knepley wrote:
On Fri, Oct 15, 2021 at 10:16 AM Pierre Seize <[email protected] <mailto:[email protected]>> 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
    <[email protected] <mailto:[email protected]>> 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
        <[email protected] <mailto:[email protected]>> 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/>

Reply via email to