fenicsxconcrete.boundary_conditions package
Submodules
fenicsxconcrete.boundary_conditions.bcs module
Easy definition of Dirichlet and Neumann BCs.
- class BoundaryConditions(domain: Mesh, space: FunctionSpace, facet_tags: ndarray | None = None)[source]
Bases:
LogMixinHandles Dirichlet and Neumann boundary conditions.
- domain
The computational domain.
- V
The finite element space.
- add_dirichlet_bc(value: Function | Constant | DirichletBC | ndarray | Callable, boundary: int | ndarray | Callable | None = None, sub: int = None, method: str = 'topological', entity_dim: int | None = None) None[source]
Adds a Dirichlet bc.
- Parameters:
value – Anything that might be used to define the Dirichlet function. It can be a Function, a Callable which is then interpolated or an already existing Dirichlet BC, or … (see type hint).
boundary – The part of the boundary whose dofs should be constrained. This can be a callable defining the boundary geometrically or an array of entity tags or an integer marking the boundary if facet_tags is not None.
sub – If sub is not None the subspace V.sub(sub) will be constrained.
method – A hint which method should be used to locate the dofs. Choices: ‘topological’ or ‘geometrical’.
entity_dim – The dimension of the entities to be located topologically. Note that entity_dim is required if sub is not None and method=geometrical.
- add_neumann_bc(marker: int, value: Constant) None[source]
Adds a Neumann BC.
- Parameters:
marker – The id of the boundary where Neumann BC should be applied.
value – The Neumann data, e.g. a traction vector. This has to be a valid ufl object.
- property bcs: list[DirichletBC]
returns the list of Dirichlet BCs
- Returns:
The list of Dirichlet BCs.
- clear(dirichlet: bool = True, neumann: bool = True) None[source]
Clears list of Dirichlet and/or Neumann BCs.
- Parameters:
dirichlet – flag for Dirichlet Bcs (if true will clear those)
neumann – flag for Neumann Bcs (if true will clear those)
- property has_dirichlet: bool
check if Dirichlet BCs are defined
- Returns:
True or False
- property has_neumann: bool
check if Neumann BCs are defined
- Returns:
True or False
- property neumann_bcs: Form
creates the ufl form of (sum of) Neumann BCs
- Returns:
A ufl object representing Neumann BCs
fenicsxconcrete.boundary_conditions.boundary module
Easy definition of boundaries.
- create_facet_tags(mesh: Mesh, boundaries: dict[str, tuple[int, Callable]]) tuple[ndarray, dict[str, int]][source]
Creates facet tags for the given mesh and boundaries.
This code is part of the FEniCSx tutorial by Jørgen S. Dokken. See https://jsdokken.com/dolfinx-tutorial/chapter3/robin_neumann_dirichlet.html?highlight=sorted_facets#implementation # noqa: E501
- Parameters:
mesh – The computational domain.
boundaries – The definition of boundaries where each key is a string and each value is a tuple of an integer and a marker function.
- Returns:
A tuple (facet_tags, marked_boundary) where facet_tags is an array with dtype int and marked_boundary is a dict where each key is a string and each value is an int.
- line_at(coordinates: list[float], dims: list[str | int]) Callable[source]
Defines a line where x[dims[0]] equals coordinates[0] and x[dims[1]] equals coordinates[1].
- Parameters:
coordinates – list of values
dims – list of dimension
- Returns:
function defining the boundary
- plane_at(coordinate: float, dim: str | int) Callable[source]
Defines a plane where x[dim] equals coordinate.
- Parameters:
coordinate – value
dim – dimension
- Returns:
function defining the boundary
- point_at(coord: Iterable[int] | Iterable[float]) Callable[source]
Defines a point.
- Parameters:
coord – points coordinates
- Returns:
function defining the boundary
- to_floats(x: Iterable[int] | Iterable[float]) list[float][source]
Converts x to a 3d coordinate.
- Parameters:
x – point coordinates at least 1D
- Returns:
point described as list with x,y,z value
- within_range(start: Iterable[int] | Iterable[float], end: Iterable[int] | Iterable[float], tol: float = 1e-06) Callable[source]
Defines a range.
It is best used together with dolfinx.mesh.locate_entities_boundary and topological definition of the Dirichlet BC, because the Callable will mark the whole range and not just the boundary.
- Parameters:
start – The start point of the range.
end – The end point of the range.
- Returns:
function defining the boundary