dolfinx.mesh
Creation, refining and marking of meshes
Functions
|
Create a :class:dolfinx.mesh.MeshTags` object that associates data with a subset of mesh entities, where the entities are defined by their vertices. |
|
Compute mesh entities satisfying a geometric marking function. |
|
Compute mesh entities that are connected to an owned boundary facet and satisfy a geometric marking function. |
|
Refine a mesh. |
|
Create a mesh from topology and geometry arrays. |
|
Create a mesh with specified entities from another mesh. |
|
Create a MeshTags object that associates data with a subset of mesh entities. |
|
|
|
|
|
Create an interval mesh. |
|
Create a mesh on the unit interval. |
|
Create a rectangle mesh. |
|
Create a mesh of a unit square. |
|
Create a box mesh. |
|
Create a mesh of a unit cube. |
|
Refine a (topologically) one dimensional mesh. |
|
Refine a mesh. |
|
Generate cell mesh tags on a refined mesh from the mesh tags on the coarse parent mesh. |
|
Compute the geometric DOFs associated with the closure of the given mesh entities. |
Classes
|
A mesh. |
|
Mesh tags associate data (markers) with a subset of mesh entities of a given dimension. |
|
|
|
- class dolfinx.mesh.CellType(value, names=<not given>, *values, module=None, qualname=None, type=None, start=1, boundary=None)
Bases:
Enum
- hexahedron = -8
- interval = 2
- property name
(self) -> object
- point = 1
- prism = -6
- pyramid = -5
- quadrilateral = -4
- tetrahedron = 4
- triangle = 3
- class dolfinx.mesh.GhostMode(value, names=<not given>, *values, module=None, qualname=None, type=None, start=1, boundary=None)
Bases:
Enum
- none = 0
- class dolfinx.mesh.Mesh(mesh, domain: Mesh | None)[source]
Bases:
object
A mesh.
Initialize mesh from a C++ mesh.
- Parameters:
mesh – A C++ mesh object.
domain – A UFL domain.
Note
Mesh objects should not usually be created using this initializer directly.
- property comm
- property geometry
Mesh geometry.
- h(dim: int, entities: ndarray[Any, dtype[int32]]) ndarray[Any, dtype[float64]] [source]
Geometric size measure of cell entities.
- Parameters:
dim – Topological dimension of the entities to compute the size measure of.
entities – Indices of entities of dimension
dim
to compute size measure of.
- Returns:
Size measure for each requested entity.
- property name
- property topology
Mesh topology.
- class dolfinx.mesh.MeshTags(meshtags)[source]
Bases:
object
Mesh tags associate data (markers) with a subset of mesh entities of a given dimension.
Initialize tags from a C++ MeshTags object.
- Parameters:
meshtags – C++ mesh tags object.
Note
MeshTags objects should not usually be created using this initializer directly.
A Python mesh is passed to the initializer as it may have UFL data attached that is not attached the C + + Mesh that is associated with the C + +
meshtags
object. If mesh is passed,mesh
andmeshtags
must share the same C + + mesh.- property dim: int
Topological dimension of the tagged entities.
- find(value) ndarray[Any, dtype[int32]] [source]
Get a list of all entity indices with a given value.
- Parameters:
value – Tag value to search for.
- Returns:
Indices of entities with tag
value
.
- property indices: ndarray[Any, dtype[int32]]
Indices of tagged mesh entities.
- property name: str
Name of the mesh tags object.
- property values
Values associated with tagged mesh entities.
- dolfinx.mesh.compute_incident_entities(topology, entities: ndarray[Any, dtype[int32]], d0: int, d1: int)[source]
- dolfinx.mesh.create_box(comm: ~mpi4py.MPI.Comm, points: list[~collections.abc.Buffer | ~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]] | ~numpy._typing._nested_sequence._NestedSequence[~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]]] | bool | int | float | complex | str | bytes | ~numpy._typing._nested_sequence._NestedSequence[bool | int | float | complex | str | bytes]], n: list, cell_type=CellType.tetrahedron, dtype: ~numpy.dtype[~typing.Any] | None | type[~typing.Any] | ~numpy._typing._dtype_like._SupportsDType[~numpy.dtype[~typing.Any]] | str | tuple[~typing.Any, int] | tuple[~typing.Any, ~typing.SupportsIndex | ~collections.abc.Sequence[~typing.SupportsIndex]] | list[~typing.Any] | ~numpy._typing._dtype_like._DTypeDict | tuple[~typing.Any, ~typing.Any] = <class 'numpy.float64'>, ghost_mode=GhostMode.shared_facet, partitioner=None) Mesh [source]
Create a box mesh.
- Parameters:
comm – MPI communicator.
points – Coordinates of the ‘lower-left’ and ‘upper-right’ corners of the box.
n – List of cells in each direction
cell_type – The cell type.
dtype – Float type for the mesh geometry(
numpy.float32
ornumpy.float64
).ghost_mode – The ghost mode used in the mesh partitioning.
partitioner – Function that computes the parallel distribution of cells across MPI ranks.
- Returns:
A mesh of a box domain.
- dolfinx.mesh.create_interval(comm: ~mpi4py.MPI.Comm, nx: int, points: ~collections.abc.Buffer | ~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]] | ~numpy._typing._nested_sequence._NestedSequence[~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]]] | bool | int | float | complex | str | bytes | ~numpy._typing._nested_sequence._NestedSequence[bool | int | float | complex | str | bytes], dtype: ~numpy.dtype[~typing.Any] | None | type[~typing.Any] | ~numpy._typing._dtype_like._SupportsDType[~numpy.dtype[~typing.Any]] | str | tuple[~typing.Any, int] | tuple[~typing.Any, ~typing.SupportsIndex | ~collections.abc.Sequence[~typing.SupportsIndex]] | list[~typing.Any] | ~numpy._typing._dtype_like._DTypeDict | tuple[~typing.Any, ~typing.Any] = <class 'numpy.float64'>, ghost_mode=GhostMode.shared_facet, partitioner=None) Mesh [source]
Create an interval mesh.
- Parameters:
comm – MPI communicator.
nx – Number of cells.
points – Coordinates of the end points.
dtype – Float type for the mesh geometry(
numpy.float32
ornumpy.float64
).ghost_mode – Ghost mode used in the mesh partitioning. Options are
GhostMode.none
andGhostMode.shared_facet
.partitioner – Partitioning function to use for determining the parallel distribution of cells across MPI ranks.
- Returns:
An interval mesh.
- dolfinx.mesh.create_mesh(comm: Comm, cells: ndarray[Any, dtype[int64]], x: ndarray[Any, dtype[floating]], e: Mesh | FiniteElement | _BasixElement | CoordinateElement, partitioner: Callable | None = None) Mesh [source]
Create a mesh from topology and geometry arrays.
- Parameters:
comm – MPI communicator to define the mesh on.
cells – Cells of the mesh.
cells[i]
are the ‘nodes’ of celli
.x – Mesh geometry (‘node’ coordinates), with shape
(num_nodes, gdim)
.e – UFL mesh. The mesh scalar type is determined by the scalar type of
e
.partitioner – Function that computes the parallel distribution of cells across MPI ranks.
Note
If required, the coordinates
x
will be cast to the same type as the domain/elemente
.- Returns:
A mesh.
- dolfinx.mesh.create_rectangle(comm: ~mpi4py.MPI.Comm, points: ~collections.abc.Buffer | ~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]] | ~numpy._typing._nested_sequence._NestedSequence[~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]]] | bool | int | float | complex | str | bytes | ~numpy._typing._nested_sequence._NestedSequence[bool | int | float | complex | str | bytes], n: ~collections.abc.Buffer | ~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]] | ~numpy._typing._nested_sequence._NestedSequence[~numpy._typing._array_like._SupportsArray[~numpy.dtype[~typing.Any]]] | bool | int | float | complex | str | bytes | ~numpy._typing._nested_sequence._NestedSequence[bool | int | float | complex | str | bytes], cell_type=CellType.triangle, dtype: ~numpy.dtype[~typing.Any] | None | type[~typing.Any] | ~numpy._typing._dtype_like._SupportsDType[~numpy.dtype[~typing.Any]] | str | tuple[~typing.Any, int] | tuple[~typing.Any, ~typing.SupportsIndex | ~collections.abc.Sequence[~typing.SupportsIndex]] | list[~typing.Any] | ~numpy._typing._dtype_like._DTypeDict | tuple[~typing.Any, ~typing.Any] = <class 'numpy.float64'>, ghost_mode=GhostMode.shared_facet, partitioner=None, diagonal: ~dolfinx.cpp.mesh.DiagonalType = DiagonalType.right) Mesh [source]
Create a rectangle mesh.
- Parameters:
comm – MPI communicator.
points – Coordinates of the lower - left and upper - right corners of the rectangle.
n – Number of cells in each direction.
cell_type – Mesh cell type.
dtype – Float type for the mesh geometry(
numpy.float32
ornumpy.float64
)ghost_mode – Ghost mode used in the mesh partitioning.
partitioner – Function that computes the parallel distribution of cells across MPI ranks.
diagonal – Direction of diagonal of triangular meshes. The options are
left
,right
,crossed
,left / right
,right / left
.
- Returns:
A mesh of a rectangle.
- dolfinx.mesh.create_submesh(msh: Mesh, dim: int, entities: ndarray[Any, dtype[int32]]) tuple[Mesh, ndarray[Any, dtype[int32]], ndarray[Any, dtype[int32]], ndarray[Any, dtype[int32]]] [source]
Create a mesh with specified entities from another mesh.
- Parameters:
mesh – Mesh to create the sub-mesh from.
dim – Topological dimension of the entities in
msh
to include in the sub-mesh.entities – Indices of entities in
msh
to include in the sub-mesh.
- Returns:
The (1) sub mesh, (2) entity map, (3) vertex map, and (4) node map (geometry). Each of the maps a local index of the sub mesh to a local index of
msh
.
- dolfinx.mesh.create_unit_cube(comm: ~mpi4py.MPI.Comm, nx: int, ny: int, nz: int, cell_type=CellType.tetrahedron, dtype: ~numpy.dtype[~typing.Any] | None | type[~typing.Any] | ~numpy._typing._dtype_like._SupportsDType[~numpy.dtype[~typing.Any]] | str | tuple[~typing.Any, int] | tuple[~typing.Any, ~typing.SupportsIndex | ~collections.abc.Sequence[~typing.SupportsIndex]] | list[~typing.Any] | ~numpy._typing._dtype_like._DTypeDict | tuple[~typing.Any, ~typing.Any] = <class 'numpy.float64'>, ghost_mode=GhostMode.shared_facet, partitioner=None) Mesh [source]
Create a mesh of a unit cube.
- Parameters:
comm – MPI communicator.
nx – Number of cells in “x” direction.
ny – Number of cells in “y” direction.
nz – Number of cells in “z” direction.
cell_type – Mesh cell type
dtype – Float type for the mesh geometry(
numpy.float32
ornumpy.float64
).ghost_mode – Ghost mode used in the mesh partitioning.
partitioner – Function that computes the parallel distribution of cells across MPI ranks.
- Returns:
- A mesh of an axis-aligned unit cube with corners at
(0, 0, 0)
and
(1, 1, 1)
.
- A mesh of an axis-aligned unit cube with corners at
- dolfinx.mesh.create_unit_interval(comm: ~mpi4py.MPI.Comm, nx: int, dtype: ~numpy.dtype[~typing.Any] | None | type[~typing.Any] | ~numpy._typing._dtype_like._SupportsDType[~numpy.dtype[~typing.Any]] | str | tuple[~typing.Any, int] | tuple[~typing.Any, ~typing.SupportsIndex | ~collections.abc.Sequence[~typing.SupportsIndex]] | list[~typing.Any] | ~numpy._typing._dtype_like._DTypeDict | tuple[~typing.Any, ~typing.Any] = <class 'numpy.float64'>, ghost_mode=GhostMode.shared_facet, partitioner=None) Mesh [source]
Create a mesh on the unit interval.
- Parameters:
comm – MPI communicator.
nx – Number of cells.
points – Coordinates of the end points.
dtype – Float type for the mesh geometry(
numpy.float32
ornumpy.float64
).ghost_mode – Ghost mode used in the mesh partitioning. Options are
GhostMode.none
andGhostMode.shared_facet
.partitioner – Partitioning function to use for determining the parallel distribution of cells across MPI ranks.
- Returns:
A unit interval mesh with end points at 0 and 1.
- dolfinx.mesh.create_unit_square(comm: ~mpi4py.MPI.Comm, nx: int, ny: int, cell_type=CellType.triangle, dtype: ~numpy.dtype[~typing.Any] | None | type[~typing.Any] | ~numpy._typing._dtype_like._SupportsDType[~numpy.dtype[~typing.Any]] | str | tuple[~typing.Any, int] | tuple[~typing.Any, ~typing.SupportsIndex | ~collections.abc.Sequence[~typing.SupportsIndex]] | list[~typing.Any] | ~numpy._typing._dtype_like._DTypeDict | tuple[~typing.Any, ~typing.Any] = <class 'numpy.float64'>, ghost_mode=GhostMode.shared_facet, partitioner=None, diagonal: ~dolfinx.cpp.mesh.DiagonalType = DiagonalType.right) Mesh [source]
Create a mesh of a unit square.
- Parameters:
comm – MPI communicator.
nx – Number of cells in the “x” direction.
ny – Number of cells in the “y” direction.
cell_type – Mesh cell type.
dtype – Float type for the mesh geometry(
numpy.float32
ornumpy.float64
).ghost_mode – Ghost mode used in the mesh partitioning.
partitioner – Function that computes the parallel distribution of cells across MPI ranks.
diagonal – Direction of diagonal.
- Returns:
A mesh of a square with corners at (0, 0) and (1, 1).
- dolfinx.mesh.entities_to_geometry(mesh: Mesh, dim: int, entities: ndarray[Any, dtype[int32]], permute=False) ndarray[Any, dtype[int32]] [source]
Compute the geometric DOFs associated with the closure of the given mesh entities.
- Parameters:
mesh – The mesh.
dim – Topological dimension of the entities of interest.
entities – Entity indices (local to the process).
permute – Permute the DOFs such that they are consistent with the orientation of dim-dimensional mesh entities. This requires create_entity_permutations to be called first.
- Returns:
The geometric DOFs associated with the closure of the entities in entities.
- dolfinx.mesh.locate_entities(mesh: Mesh, dim: int, marker: Callable) ndarray [source]
Compute mesh entities satisfying a geometric marking function.
- Parameters:
mesh – Mesh to locate entities on.
dim – Topological dimension of the mesh entities to consider.
marker – A function that takes an array of points
x
with shape(gdim, num_points)
and returns an array of booleans of lengthnum_points
, evaluating to True for entities to be located.
- Returns:
Indices (local to the process) of marked mesh entities.
- dolfinx.mesh.locate_entities_boundary(mesh: Mesh, dim: int, marker: Callable) ndarray [source]
Compute mesh entities that are connected to an owned boundary facet and satisfy a geometric marking function.
For vertices and edges, in parallel this function will not necessarily mark all entities that are on the exterior boundary. For example, it is possible for a process to have a vertex that lies on the boundary without any of the attached facets being a boundary facet. When used to find degrees-of-freedom, e.g. using
dolfinx.fem.locate_dofs_topological()
, the function that uses the data returned by this function must typically perform some parallel communication.- Parameters:
mesh – Mesh to locate boundary entities on.
dim – Topological dimension of the mesh entities to consider
marker – Function that takes an array of points
x
with shape(gdim, num_points)
and returns an array of booleans of lengthnum_points
, evaluating toTrue
for entities to be located.
- Returns:
Indices (local to the process) of marked mesh entities.
- dolfinx.mesh.meshtags(mesh: Mesh, dim: int, entities: ndarray[Any, dtype[int32]], values: ndarray | int | float) MeshTags [source]
Create a MeshTags object that associates data with a subset of mesh entities.
- Parameters:
mesh – The mesh.
dim – Topological dimension of the mesh entity.
entities – Indices(local to process) of entities to associate values with . The array must be sorted and must not contain duplicates.
values – The corresponding value for each entity.
- Returns:
A mesh tags object.
Note
The type of the returned MeshTags is inferred from the type of
values
.
- dolfinx.mesh.meshtags_from_entities(mesh: Mesh, dim: int, entities: AdjacencyList_int32, values: ndarray[Any, dtype[Any]])[source]
Create a :class:dolfinx.mesh.MeshTags` object that associates data with a subset of mesh entities, where the entities are defined by their vertices.
- Parameters:
mesh – The mesh.
dim – Topological dimension of the mesh entity.
entities – Entities to associated values with, with entities defined by their vertices.
values – The corresponding value for each entity.
- Returns:
A mesh tags object.
Note
The type of the returned MeshTags is inferred from the type of
values
.
- dolfinx.mesh.refine(mesh: Mesh, edges: ndarray | None = None, redistribute: bool = True) Mesh [source]
Refine a mesh.
- Parameters:
mesh – Mesh from which to create the refined mesh.
edges – Indices of edges to split during refinement. If
None
, mesh refinement is uniform.redistribute – Refined mesh is re-partitioned if
True
.
- Returns:
Refined mesh.
- dolfinx.mesh.refine_interval(mesh: Mesh, cells: ndarray | None = None, redistribute: bool = True, ghost_mode: GhostMode = GhostMode.shared_facet) tuple[Mesh, ndarray[Any, dtype[int32]]] [source]
Refine a (topologically) one dimensional mesh.
- Parameters:
mesh – Mesh to refine
cells – Indices of cells, i.e. edges, to split druing refinement. If
None
, mesh refinement is uniform.redistribute – Refined mesh is re-partitioned if
True
.ghost_mode – ghost mode of the refined mesh
- Returns:
Refined mesh and parent cells
- dolfinx.mesh.refine_plaza(mesh: Mesh, edges: ndarray | None = None, redistribute: bool = True, option: RefinementOption = RefinementOption.none) tuple[Mesh, ndarray[Any, dtype[int32]], ndarray[Any, dtype[int32]]] [source]
Refine a mesh.
- Parameters:
mesh – Mesh from which to create the refined mesh.
edges – Indices of edges to split during refinement. If
None
, mesh refinement is uniform.redistribute – Refined mesh is re-partitioned if
True
.option – Control computation of the parent-refined mesh data.
- Returns:
Refined mesh, list of parent cell for each refine cell, and list
- dolfinx.mesh.transfer_meshtag(meshtag: MeshTags, mesh1: Mesh, parent_cell: ndarray[Any, dtype[int32]], parent_facet: ndarray[Any, dtype[int8]] | None = None) MeshTags [source]
Generate cell mesh tags on a refined mesh from the mesh tags on the coarse parent mesh.
- Parameters:
meshtag – Mesh tags on the coarse, parent mesh.
mesh1 – The refined mesh.
parent_cell – Index of the parent cell for each cell in the refined mesh.
parent_facet – Index of the local parent facet for each cell in the refined mesh. Only required for transfer tags on facets.
- Returns:
Mesh tags on the refined mesh.