petsc-3.11.4 2019-09-28
Report Typos and Errors

DMPlexCreateFromDAG

This takes as input the adjacency-list representation of the Directed Acyclic Graph (Hasse Diagram) encoding a mesh, and produces a DM

Synopsis

#include "petscdmplex.h"   
#include "petscdmplex.h"   
PetscErrorCode DMPlexCreateFromDAG(DM dm, PetscInt depth, const PetscInt numPoints[], const PetscInt coneSize[], const PetscInt cones[], const PetscInt coneOrientations[], const PetscScalar vertexCoords[])

Input Parameters

dm - The empty DM object, usually from DMCreate() and DMSetDimension()
depth - The depth of the DAG
numPoints - Array of size depth + 1 containing the number of points at each depth
coneSize - The cone size of each point
cones - The concatenation of the cone points for each point, the cone list must be oriented correctly for each point
coneOrientations - The orientation of each cone point
vertexCoords - An array of numPoints[0]*spacedim numbers representing the coordinates of each vertex, with spacedim the value set via DMSetCoordinateDim()

Output Parameter

dm -The DM

Note: Two triangles sharing a face would have input

 depth = 1, numPoints = [4 2], coneSize = [3 3 0 0 0 0]
 cones = [2 3 4  3 5 4], coneOrientations = [0 0 0  0 0 0]
vertexCoords = [-1.0 0.0  0.0 -1.0  0.0 1.0  1.0 0.0]

which would result in the DMPlex

       4
     / | \
    /  |  \
   /   |   \
  2  0 | 1  5
   \   |   /
    \  |  /
     \ | /
       3

Notice that all points are numbered consecutively, unlikely DMPlexCreateFromCellList()

See Also

DMPlexCreateFromCellList(), DMPlexCreate()

Level

advanced

Location

src/dm/impls/plex/plexcreate.c
Index of all DMPLEX routines
Table of Contents for all manual pages
Index of all manual pages