Mesh

class dolfin.cpp.mesh.Mesh(*args, **kwargs)

Bases: dolfin.cpp.common.Variable, dolfin.cpp.mesh.HierarchicalMesh

A Mesh consists of a set of connected and numbered mesh entities.

Both the representation and the interface are dimension-independent, but a concrete interface is also provided for standard named mesh entities:

Entity Dimension Codimension
Vertex 0  
Edge 1  
Face 2  
Facet   1
Cell   0

When working with mesh iterators, all entities and connectivity are precomputed automatically the first time an iterator is created over any given topological dimension or connectivity.

Note that for efficiency, only entities of dimension zero (vertices) and entities of the maximal dimension (cells) exist when creating a Mesh. Other entities must be explicitly created by calling init(). For example, all edges in a mesh may be created by a call to mesh.init(1). Similarly, connectivities such as all edges connected to a given vertex must also be explicitly created (in this case by a call to mesh.init(0, 1)).

Create a mesh from a filename or a geometry.

A mesh may be created from a given filename, which should contain mesh data stored in DOLFIN XML format:

mesh = Mesh("mesh.xml")

An empty mesh may be created as follows:

mesh = Mesh()

A copy of a mesh may be created as follows:

mesh_copy = Mesh(mesh)

A mesh may also be created from a CSG geometry and a resolution:

mesh = Mesh(geometry, resolution)
bounding_box_tree()

Get bounding box tree for mesh. The bounding box tree is initialized and built upon the first call to this function. The bounding box tree can be used to compute collisions between the mesh and other objects. It is the responsibility of the caller to use (and possibly rebuild) the tree. It is stored as a (mutable) member of the mesh to enable sharing of the bounding box tree data structure.

cell_orientations()

Get the cell orientations set.

Returns
numpy.array(int)
Cell orientations
cells()

Get cell connectivity.

Returns
numpy.array(int)
Connectivity for all cells.
Example
>>> mesh = dolfin.UnitSquare(1,1)
>>> mesh.cells()
array([[0, 1, 3],
      [0, 2, 3]])
clean()

Clean out all auxiliary topology data. This clears all topological data, except the connectivity between cells and vertices.

clear()

Clear all mesh data.

color()

Overloaded versions

  • color(coloring_type)

    Color the cells of the mesh such that no two neighboring cells share the same color. A colored mesh keeps a CellFunction<std::size_t> named “cell colors” as mesh data which holds the colors of the mesh.

    Arguments
    coloring_type (str)

    Coloring type, specifying what relation makes two cells neighbors, can be one of “vertex”, “edge” or “facet”.

    Returns
    numpy.array(int)

    The colors as a mesh function over the cells of the mesh.

  • color(coloring_type)

    Color the cells of the mesh such that no two neighboring cells share the same color. A colored mesh keeps a CellFunction<std::size_t> named “cell colors” as mesh data which holds the colors of the mesh.

    Arguments
    coloring_type (numpy.array(int))

    Coloring type given as list of topological dimensions, specifying what relation makes two mesh entinties neighbors.

    Returns
    numpy.array(int)

    The colors as a mesh function over entities of the mesh.

coordinates()
  • coordinates()

    Get vertex coordinates.

    Returns
    numpy.array(float)

    Coordinates of all vertices.

    Example
    >>> mesh = dolfin.UnitSquare(1,1)
    >>> mesh.coordinates()
    array([[ 0.,  0.],
           [ 1.,  0.],
           [ 0.,  1.],
           [ 1.,  1.]])
    
data()

Overloaded versions

  • data()

    Get mesh data.

    Returns
    MeshData

    The mesh data object associated with the mesh.

  • data()

    Get mesh data (const version).

domains()

Overloaded versions

  • domains()

    Get mesh (sub)domains.

    Returns
    MeshDomains

    The (sub)domains associated with the mesh.

  • domains()

    Get mesh (sub)domains.

geometry()

Overloaded versions

  • geometry()

    Get mesh geometry.

    Returns
    MeshGeometry

    The geometry object associated with the mesh.

  • geometry()

    Get mesh geometry (const version).

hash()

Compute hash of mesh, currently based on the has of the mesh geometry and mesh topology.

Returns
int
A tree-hashed value of the coordinates over all MPI processes
hmax()

Compute maximum cell diameter.

Returns
float
The maximum cell diameter, the diameter is computed as two times the circumradius (http://mathworld.wolfram.com).
Example
>>> mesh = dolfin.UnitSquare(2,2)
>>> mesh.hmax()
0.70710678118654757
hmin()

Compute minimum cell diameter.

Returns
float
The minimum cell diameter, the diameter is computed as two times the circumradius (http://mathworld.wolfram.com).
Example
>>> mesh = dolfin.UnitSquare(2,2)
>>> mesh.hmin()
0.70710678118654757
init()

Overloaded versions

  • init(dim)

    Compute entities of given topological dimension.

    Arguments
    dim (int)

    Topological dimension.

    Returns
    int

    Number of created entities.

  • init(d0, d1)

    Compute connectivity between given pair of dimensions.

    Arguments
    d0 (int)

    Topological dimension.

    d1 (int)

    Topological dimension.

  • init()

    Compute all entities and connectivity.

init_cell_orientations()

Compute and initialize cell_orientations relative to a given global outward direction/normal/orientation. Only defined if mesh is orientable.

Arguments
global_normal (Expression)
A global normal direction to the mesh
move()

Overloaded versions

  • move(boundary)

    Move coordinates of mesh according to new boundary coordinates.

    Arguments
    boundary (BoundaryMesh)

    A mesh containing just the boundary cells.

    Returns
    MeshDisplacement

    Displacement encapsulated in Expression subclass MeshDisplacement.

  • move(mesh)

    Move coordinates of mesh according to adjacent mesh with common global vertices.

    Arguments
    mesh (Mesh)

    A Mesh object.

    Returns
    MeshDisplacement

    Displacement encapsulated in Expression subclass MeshDisplacement.

  • move(displacement)

    Move coordinates of mesh according to displacement function.

    Arguments
    displacement (GenericFunction)

    A GenericFunction object.

mpi_comm()

Mesh MPI communicator

num_cells()

Get number of cells in mesh.

Returns
int
Number of cells.
Example

Note

No example code available for this function.

num_edges()

Get number of edges in mesh.

Returns
int
Number of edges.
Example

Note

No example code available for this function.

num_entities()

Get number of entities of given topological dimension.

Arguments
d (int)
Topological dimension.
Returns
int
Number of entities of topological dimension d.
Example

Note

No example code available for this function.

num_faces()

Get number of faces in mesh.

Returns
int
Number of faces.
Example

Note

No example code available for this function.

num_facets()

Get number of facets in mesh.

Returns
int
Number of facets.
Example

Note

No example code available for this function.

num_vertices()

Get number of vertices in mesh.

Returns
int
Number of vertices.
Example

Note

No example code available for this function.

order()

Order all mesh entities.

See also

UFC documentation (put link here!)

ordered()

Check if mesh is ordered according to the UFC numbering convention.

Returns
bool
The return values is true iff the mesh is ordered.
renumber_by_color()
rmax()

Compute maximum cell inradius.

Returns
float
The maximum of cells’ inscribed sphere radii
Example

Note

No example code available for this function.

rmin()

Compute minimum cell inradius.

Returns
float
The minimum of cells’ inscribed sphere radii
Example

Note

No example code available for this function.

rotate()

Overloaded versions

  • rotate(angle, axis=2)

    Rotate mesh around a coordinate axis through center of mass of all mesh vertices

    Arguments
    angle (float)

    The number of degrees (0-360) of rotation.

    axis (int)

    The coordinate axis around which to rotate the mesh.

  • rotate(angle, axis, point)

    Rotate mesh around a coordinate axis through a given point

    Arguments
    angle (float)

    The number of degrees (0-360) of rotation.

    axis (int)

    The coordinate axis around which to rotate the mesh.

    point (Point)

    The point around which to rotate the mesh.

size()

Get number of local entities of given topological dimension.

Arguments
dim (int)
Topological dimension.
Returns
int
Number of local entities of topological dimension d.
Example

Note

No example code available for this function.

size_global()

Get global number of entities of given topological dimension.

Arguments
dim (int)
Topological dimension.
Returns
int
Global number of entities of topological dimension d.
Example

Note

No example code available for this function.

smooth()

Smooth internal vertices of mesh by local averaging.

Arguments
num_iterations (int)
Number of iterations to perform smoothing, default value is 1.
smooth_boundary()

Smooth boundary vertices of mesh by local averaging.

Arguments
num_iterations (int)
Number of iterations to perform smoothing, default value is 1.
harmonic_smoothing (bool)
Flag to turn on harmonics smoothing, default value is true.
snap_boundary()

Snap boundary vertices of mesh to match given sub domain.

Arguments
sub_domain (SubDomain)
A SubDomain object.
harmonic_smoothing (bool)
Flag to turn on harmonics smoothing, default value is true.
thisown

The membership flag

topology()

Overloaded versions

  • topology()

    Get mesh topology.

    Returns
    MeshTopology

    The topology object associated with the mesh.

  • topology()

    Get mesh topology (const version).

translate()

Translate mesh according to a given vector.

Arguments
point (Point)
The vector defining the translation.
type()

Overloaded versions

  • type()

    Get mesh cell type.

    Returns
    CellType

    The cell type object associated with the mesh.

  • type()

    Get mesh cell type (const version).

ufl_cell()

Returns the ufl cell of the mesh.

The cell corresponds to the topological dimension of the mesh.

ufl_domain()

Returns the ufl Domain corresponding to the mesh.

ufl_id()

Returns an id that UFL can use to decide if two objects are the same.