API reference#
- class scifem.PointSource(V: FunctionSpace, points: ndarray[tuple[int, ...], dtype[float32]] | ndarray[tuple[int, ...], dtype[float64]], magnitude: floating | complexfloating = np.float64(1.0))[source]#
Class for defining a point source in a given function space.
- apply_to_vector(b: Function | Vector | Vec, recompute: bool = False)[source]#
Apply the point sources to a vector.
- Parameters:
b – The vector to apply the point sources to.
recompute – If the point sources should be recomputed before applying. Recomputation should only be done if the mesh geometry has been modified.
Note
The user is responsible for forward scattering of the vector after applying the point sources.
Note
If a PETSc vector is passed in, one has to call
b.assemble()
prior to solving the linear system (post scattering).
- scifem.assemble_scalar(J: Form | Form, entity_maps: dict[Mesh, ndarray[tuple[int, ...], dtype[int32]]] | None = None) floating | complexfloating [source]#
Assemble a scalar form and gather result across processes
- Parameters:
form – The form to assemble.
entity_maps – Maps of entities on related submeshes to the domain used in J.
- Returns:
The accumulated value of the assembled form.
- scifem.compute_interface_data(cell_tags: MeshTags, facet_indices: ndarray[tuple[int, ...], dtype[int32]]) ndarray[tuple[int, ...], dtype[int32]] [source]#
Compute interior facet integrals that are consistently ordered according to the cell_tags, such that the data (cell0, facet_idx0, cell1, facet_idx1) is ordered such that cell_tags[cell0]`<`cell_tags[cell1], i.e the cell with the lowest cell marker is considered the “+” restriction”.
- Parameters:
cell_tags – MeshTags that must contain an integer marker for all cells adjacent to the facet_indices
facet_indices – List of facets (local index) that are on the interface.
- Returns:
The integration data.
- scifem.compute_subdomain_exterior_facets(mesh: Mesh, ct: MeshTags, markers: Sequence[int]) ndarray[tuple[int, ...], dtype[int32]] [source]#
Find the the facets that are considered to be on the “exterior” boundary of a subdomain.
The subdomain is defined as the collection of cells in
ct
that is marked with any of themarkers
. The exterior boundary of the subdomain is defined as the collection of facets that are only connected to a single cell within the subdomain.Note
Ghosted facets are included in the resulting array.
- Parameters:
mesh – Mesh to extract subdomains from
ct – MeshTags object marking subdomains
markers – The tags making up the “new” mesh
- Returns:
The exterior facets
- scifem.create_entity_markers(domain: dolfinx.mesh.Mesh, dim: int, entities_list: list[TaggedEntities]) dolfinx.mesh.MeshTags [source]#
Mark entities of specified dimension according to a geometrical marker function.
- Parameters:
domain – A
dolfinx.mesh.Mesh
objectdim – Dimension of the entities to mark
entities_list –
A list of tuples with the following elements:
index 0
: The tag to assign to the entitiesindex 1
: A function that takes a point and returns a boolean array indicating whether the point is inside the entityindex 2
: Optional, if True, the entities will be marked on the boundary
- Returns:
A
dolfinx.mesh.MeshTags
object with the corresponding entities marked. If an entity satisfies multiple input marker functions, it is not deterministic what value the entity gets.
Note:
- scifem.create_real_functionspace(mesh: Mesh, value_shape: tuple[int, ...] = ()) FunctionSpace [source]#
Create a real function space.
- Parameters:
mesh – The mesh the real space is defined on.
value_shape – The shape of the values in the real space.
- Returns:
The real valued function space.
Note
For scalar elements value shape is
()
.
- scifem.create_space_of_simple_functions(mesh: Mesh, cell_tag: MeshTags, tags: Sequence[int | int32 | int64], value_shape: tuple[()] | tuple[int | tuple[int]] | None = None) FunctionSpace [source]#
Create a space of simple functions.
This is a space that represents piecewise constant functions of N patches, where N is the number of input tags. Each patch is defined by the cells in cell_tag which is marked with the corresponding tag value.
Note
All cells are expected to have a tag in tags.
- Parameters:
mesh – The mesh the function space is defined on.
cell_tag – The mesh tags defining the different patches of cells.
tags – The set of unique values within cell_tag that defines the different patches.
value_shape – The shape of the values in the space.
- Returns:
The space of simple functions.
- scifem.dof_to_vertexmap(V: FunctionSpace) ndarray[tuple[int, ...], dtype[int32]] [source]#
Create a map from the degrees of freedom to the vertices of the mesh. As not every degree of freedom is associated with a vertex, every dof that is not associated with a vertex returns -1
- Parameters:
V – The function space
- Returns:
An array mapping local dof i to a local vertex
- scifem.evaluate_function(u: Function, points: Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | bool | int | float | complex | str | bytes | _NestedSequence[bool | int | float | complex | str | bytes], broadcast=True) ndarray[tuple[int, ...], dtype[float64]] [source]#
Evaluate a function at a set of points.
- Parameters:
u – The function to evaluate.
points – The points to evaluate the function at.
broadcast –
If True, the values will be broadcasted to all processes.
Note
Uses a global MPI call to broadcast values, thus this has to be called on all active processes synchronously.
Note
If the function is discontinuous, different processes may return different values for the same point. In this case, the value returned is the maximum value across all processes.
- Returns:
The values of the function evaluated at the points.
- scifem.extract_submesh(mesh: Mesh, entity_tag: MeshTags, tags: Sequence[int]) SubmeshData [source]#
Generate a sub-mesh from a subset of tagged entities in a meshtag object.
- Parameters:
mesh – The mesh to extract the submesh from.
entity_tag – MeshTags object containing marked entities.
tags – What tags the marked entities used in the submesh should have.
- Returns:
A tuple (submesh, subcell_to_parent_entity, subvertex_to_parent_vertex, subnode_to_parent_node, entity_tag_on_submesh).
- scifem.find_interface(cell_tags: MeshTags, id_0: tuple[int, ...], id_1: tuple[int, ...]) ndarray[tuple[int, ...], dtype[int32]] [source]#
Given to sets of cells, find the facets that are shared between them.
- Parameters:
cell_tags – MeshTags object marking cells.
id_0 – Tags to extract for domain 0
id_1 – Tags to extract for domain 1
- Returns:
The facets shared between the two domains.
- scifem.interpolate_function_onto_facet_dofs(Q: FunctionSpace, expr: Expr, facets: ndarray[tuple[int, ...], dtype[int32]]) Function [source]#
Create a function \(u_h\in Q\) such that \(u_h=\text{expr}\) for all dofs belonging to a subset of
facets
. All other dofs are set to zero.Note
The resulting function is only correct in the “normal” direction, i.e. \(u_{bc}\cdot n = expr\), while the tangential component is uncontrolled. This makes it hard to visualize the function when outputting it to file, either through interpolation to an appropriate DG space, or to a point-cloud.
- Parameters:
Q – The function space to create the function $u_h$ in.
expr – The expression to evaluate.
facets – The facets on which to evaluate the expression.
- scifem.norm(expr: Expr, norm_type: Literal['L2', 'H1', 'H10'], entity_maps: dict[Mesh, ndarray[tuple[int, ...], dtype[int32]]] | None = None) Form [source]#
Compile the norm of an UFL expression into a DOLFINx form.
- Parameters:
expr – UFL expression
norm_type – Type of norm
entity_maps – Mapping for Constants and Coefficients within the expression that lives on a submesh.
- scifem.reverse_mark_entities(entity_map: IndexMap, entities: ndarray[tuple[int, ...], dtype[int32]]) ndarray[tuple[int, ...], dtype[int32]] [source]#
Communicate entities marked on a single process to all processes that ghosts or owns this entity.
- Parameters:
entity_map – Index-map describing entity ownership
entities – Local indices of entities to communicate
- Returns:
Local indices marked on any process sharing this entity
- scifem.transfer_meshtags_to_submesh(entity_tag: MeshTags, submesh: Mesh, vertex_to_parent: _EntityMap | ndarray[tuple[int, ...], dtype[int32]], cell_to_parent: _EntityMap | ndarray[tuple[int, ...], dtype[int32]]) tuple[MeshTags, ndarray[tuple[int, ...], dtype[int32]]] [source]#
Transfer a
entity_tag
from a parent mesh to asubmesh
.- Parameters:
entity_tag – Tag to transfer
submesh – Submesh to transfer tag to
vertex_to_parent – Mapping from submesh vertices to parent mesh vertices
cell_to_parent – Mapping from submesh cells to parent entities
- Returns:
submesh_tag
is the tag on the submesh andsub_to_parent_entity_map
is a mapping from submesh entities in the tag to the corresponding entities in the parent.- Return type:
A tuple (submesh_tag, sub_to_parent_entity_map) where
- scifem.vertex_to_dofmap(V: FunctionSpace) ndarray[tuple[int, ...], dtype[int32]] [source]#
Create a map from the vertices (local to the process) to the correspondning degrees of freedom.
- Parameters:
V – The function space
- Returns:
An array mapping local vertex i to local degree of freedom
Note
If using a blocked space this map is not unrolled for the DofMap block size.
- scifem.petsc.apply_lifting_and_set_bc(b: Vec, a: Iterable[Form] | Iterable[Iterable[Form]], bcs: Iterable[DirichletBC] | Iterable[Iterable[DirichletBC]], x: Vec | None = None, alpha: float = 1.0)[source]#
Apply lifting to a vector and set boundary conditions.
Convenience function to apply lifting and set boundary conditions for multiple matrix types. This modifies the vector b such that
\[\begin{split}b = \begin{pmatrix} b_{free} \\ b_{bc} \end{pmatrix} = \begin{pmatrix} b_{free} - \alpha \sum_{i=0}^n a[i] (u_{bc}[i] - x[i])\\ u_{bc}[0]\\ \vdots\\ u_{bc}[n] \end{pmatrix}.\end{split}\]where \(b_{free}\) is the free part of the vector, \(b_{bc}\) is the part that has boundary conditions applied, \(u_{bc}[i]\) is the value of the ith boundary condition.
- Parameters:
b – The vector to apply lifting to.
a – Sequence fo forms to apply lifting from. If the system is blocked or nested, his is a nested list of forms.
bcs – The boundary conditions to apply. If the form is blocked or nested, this is a list, while if it is a single form, this is a nested list.
x – Vector to subtract from the boundary conditions. Usually used in a Newton iteration.
alpha – The scaling factor for the boundary conditions.