[petsc-users] orientation of DMPLEX for tetrehedral meshes

Matteo Semplice matteo.semplice at uninsubria.it
Fri Jul 3 09:39:45 CDT 2026


Thanks both for your replies!

Matteo

Il 03/07/26 01:56, Matthew Knepley ha scritto:
> On Thu, Jul 2, 2026 at 10:55 AM Matteo Semplice via petsc-users 
> <petsc-users at mcs.anl.gov> wrote:
>
>     Dear PETSc,
>
>         I am a bit puzzled by the 3d tet mesh representations as
>     DMPlex. I boiled it down to a mesh of a single tet. More specifically,
>
>       * if a I create the mesh as
>
>               const PetscInt dim=3;
>               const PetscInt numCells=1; // The number of cells, only
>             on process 0
>               const PetscInt numVertices=4; // The number of vertices
>             owned by this process, or PETSC_DECIDE, only on process 0
>               const PetscInt numCorners=4; // The number of vertices
>             for each cell, only on process 0
>               const PetscBool interpolate=PETSC_TRUE; // should
>             intermediate mesh entities (faces, edges) be created?
>               // An array of numCells×numCornersnumCells×numCorners
>             numbers, the vertices for each cell, only on process 0
>               const PetscInt cells[] = {0,1,2,3 };
>               const PetscInt spaceDim=3; // The spatial dimension used
>             for coordinates
>               // An array of numVertices×spaceDimnumVertices×spaceDim
>             numbers, the coordinates of each vertex, only on process 0
>               const PetscReal vertexCoords[] = { 0,0,0,
>                                                  1,0,0,
>                                                  0,1,0,
>                                                  0,0,1 };
>               PetscCall( DMPlexCreateFromCellListPetsc(MPI_COMM_WORLD,
>             dim,
>                          numCells, numVertices, numCorners,
>             interpolate, cells,
>                          spaceDim, vertexCoords, dmMesh) );
>
>         I get the representation that I would expect:
>
>         DM Object: Mesh 1 MPI process
>          type: plex
>         Mesh in 3 dimensions:
>          Number of 0-cells per rank: 4
>          Number of 1-cells per rank: 6
>          Number of 2-cells per rank: 4
>          Number of 3-cells per rank: 1
>         Labels:
>          celltype: 4 strata with value/size (0 (4), 1 (6), 3 (4), 6 (1))
>          depth: 4 strata with value/size (0 (4), 1 (6), 2 (4), 3 (1))
>         Mesh with 1 cells
>
>         with DMPlexGetCellCoordinates giving me the coordinates
>
>         CELL 0 with vertices at
>         ( 0.000,     0, 0.000 )
>         ( 1.000,     0, 0.000 )
>         ( 0.000,     1, 0.000 )
>         ( 0.000,     0, 1.000 )
>
>       * on the other hand, if I load the MSH file
>
>             $MeshFormat
>             2.2 0 8
>             $EndMeshFormat
>             $Nodes
>             4
>             1 0.0 0.0 0.0
>             2 1.0 0.0 0.0
>             3 0.0 1.0 0.0
>             4 0.0 0.0 1.0
>             $EndNodes
>             $Elements
>             1
>             1 4 2 0 1 1 2 3 4
>             $EndElements
>
>         with DMPlexCreateGmshFromFile, the first and second vertices
>         are swapped and I get
>
>         DM Object: Mesh 1 MPI process
>          type: plex
>         Mesh in 3 dimensions:
>          Number of 0-cells per rank: 4
>          Number of 1-cells per rank: 6
>          Number of 2-cells per rank: 4
>          Number of 3-cells per rank: 1
>         Labels:
>          celltype: 4 strata with value/size (0 (4), 1 (6), 3 (4), 6 (1))
>          depth: 4 strata with value/size (0 (4), 1 (6), 2 (4), 3 (1))
>          Cell Sets: 1 strata with value/size (0 (1))
>         Mesh with 1 cells
>
>         with DMPlexGetCellCoordinates giving me
>
>         CELL 0 with vertices at
>         ( 1.000,     0, 0.000 )
>         ( 0.000,     0, 0.000 )
>         ( 0.000,     1, 0.000 )
>         ( 0.000,     0, 1.000 )
>
>     I understand that this behaviour is somehow consistent with the
>     formulas in DMPlexComputeTetrahedronGeometry_Internal, where you
>     compute the Jacobians swappnig the order of vertices (and with 0.5
>     factors that I do not understand at all... it looks as if your
>     reference element was not the unit simplex)
>
>         if (J) {
>              for (d = 0; d < dim; d++) {
>                /* I orient with outward face normals */
>                J[d * dim + 0] = 0.5 * (PetscRealPart(coords[2 * dim +
>         d]) - PetscRealPart(coords[0 * dim + d]));
>                J[d * dim + 1] = 0.5 * (PetscRealPart(coords[1 * dim +
>         d]) - PetscRealPart(coords[0 * dim + d]));
>                J[d * dim + 2] = 0.5 * (PetscRealPart(coords[3 * dim +
>         d]) - PetscRealPart(coords[0 * dim + d]));
>              }
>         }
>
>     but still I am worried by this fact that different ways of
>     creating a mesh gives a DMPlex with different orientations.
>
>     Can you explain the rationale behind this or refer me to some
>     documentation?
>
> 1. DMPlexCreateFromCellListPetsc() is very low-level, so I let you do 
> whatever you want, including tell me to use inverted 
> cells. DMPlexCreateGmshFromFile() is guaranteed to get a properly 
> oriented simplex from GMsh, so I guarantee that it produces a properly 
> oriented Plex simplex.
>
> 2. Here is the order of vertices, and face vertices, in a Plex 
> tetrahedron.
>
> https://urldefense.us/v3/__https://gitlab.com/petsc/petsc/-/blob/main/src/dm/impls/plex/plexinterpolate.c?ref_type=heads*L194__;Iw!!G_uCfscf7eWS!d6tvrzxI6yh_hzS8CpBk5uE1ZZZurd6twtW_RnzPfkPvxUkzTAf6W2LI-clrnlasaru3zm6h5u925ioefyQ5XBiI3f4JQtzhNZ2Fkg$ 
>
> I truly believe the only sensible way to describe a simple is to have 
> every normal outward facing. The convention in Gmsh is to have some 
> normals inward and some outward, which I think is ridiculous. It 
> breaks the inherent symmetry of the simplex.
>
> 3. The reference simplex has coordinates 1 and -1, and thus volume 
> 4/3. This is the convention of FIAT 
> (https://urldefense.us/v3/__https://www.firedrakeproject.org/fiat/index.html__;!!G_uCfscf7eWS!d6tvrzxI6yh_hzS8CpBk5uE1ZZZurd6twtW_RnzPfkPvxUkzTAf6W2LI-clrnlasaru3zm6h5u925ioefyQ5XBiI3f4JQtyQtGQqeg$ )
>
> 4. You can always get my orientation for things by calling DMPlexOrient().
>
> 5. Orientations are representors for elements of the dihedral group of 
> the k-cell. I am almost done with a description of this for every 
> geometric cell supported (of course you can define any cell you want 
> in Plex). I can show you this if you want, although it is unfinished.
>
>   Thanks,
>
>      Matt
>
>     Best
>
>         Matteo
>
>     -- 
>     Prof. Matteo Semplice
>     Università degli Studi dell’Insubria
>     Dipartimento di Scienza e Alta Tecnologia – DiSAT
>     Professore Associato
>     Via Valleggio, 11 – 22100 Como (CO) – Italia
>     tel.: +39 031 2386316
>
>
>
> -- 
> 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!d6tvrzxI6yh_hzS8CpBk5uE1ZZZurd6twtW_RnzPfkPvxUkzTAf6W2LI-clrnlasaru3zm6h5u925ioefyQ5XBiI3f4JQtx1gAf7kA$  
> <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!d6tvrzxI6yh_hzS8CpBk5uE1ZZZurd6twtW_RnzPfkPvxUkzTAf6W2LI-clrnlasaru3zm6h5u925ioefyQ5XBiI3f4JQtx2d472Bw$ >

-- 
Prof. Matteo Semplice
Università degli Studi dell’Insubria
Dipartimento di Scienza e Alta Tecnologia – DiSAT
Professore Associato
Via Valleggio, 11 – 22100 Como (CO) – Italia
tel.: +39 031 2386316
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20260703/5ed3c3e8/attachment-0001.html>


More information about the petsc-users mailing list