[petsc-users] Element connectivity of a DMPlex

[Matthew Knepley]

On Thu, May 22, 2025 at 12:25 PM Noam T. <dontbugthedevs at proton.me> wrote:

> Hello,
>
> Thank you the various options.
>
> Use case here would be obtaining the exact output generated by option 1),
> DMView() with PETSC_VIEWER_HDF5_VIZ; in particular, the matrix generated
> under /viz/topology/cells.
>
> There are several ways you might do this. It helps to know what you are
> aiming for.
>
> 1) If you just want this output, it might be easier to just DMView() with
> the PETSC_VIEWER_HDF5_VIZ format, since that just outputs the cell-vertex
> topology and coordinates
>
>
> Is it possible to get this information in memory, onto a Mat, Vec or some
> other Int array object directly? it would be handy to have it in order to
> manipulate it and/or save it to a different format/file. Saving to an HDF5
> and loading it again seems redundant.
>
>
> 2) You can call DMPlexUninterpolate() to produce a mesh with just cells
> and vertices, and output it in any format.
>
> 3) If you want it in memory, but still with global indices (I don't
> understand this use case), then you can use DMPlexCreatePointNumbering()
> for an overall global numbering, or DMPlexCreateCellNumbering() and
> DMPlexCreateVertexNumbering() for separate global numberings.
>
>
> Perhaps I missed it, but getting the connectivity matrix in
> /viz/topology/cells/ did not seem directly trivial to me from the list of
> global indices returned by  DMPlexGetCell/Point/VertexNumbering() (i.e. I
> assume all the operations done when calling DMView()).
>

Something like

DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);
DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd);
DMPlexGetVertexNumbering(dm, &globalVertexNumbers);
ISGetIndices(globalVertexNumbers, &gv);
for (PetscInt c = cStart; c < cEnd; ++c) {
    PetscInt *closure = NULL;

    DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &Ncl, &closure);
    for (PetscInt cl = 0; c < Ncl * 2; cl += 2) {
        if (closure[cl] < vStart || closure[cl] >= vEnd) continue;
        const PetscInt v = gv[closure[cl]] < 0 ? -(gv[closure[cl]] + 1) :
gv[closure[cl]];

        // Do something with v
    }
    DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &Ncl, &closure);
}
ISRestoreIndices(globalVertexNumbers, &gv);
ISDestroy(&globalVertexNumbers);

  Thanks,

     Matt



> Thanks,
> Noam.
>


-- 
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://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!c-kIsacyZtVK0QT5hSQ6TYjx-_vFWjrr6SXhV0Q8o9GcMLSLwko9mlT7CXjpNSVCT4bRPKeSkCoQJMg2SkFk$  <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!c-kIsacyZtVK0QT5hSQ6TYjx-_vFWjrr6SXhV0Q8o9GcMLSLwko9mlT7CXjpNSVCT4bRPKeSkCoQJGpFQC8u$ >
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20250522/c15245d1/attachment.html>