CCFEMFunctionSpace

class dolfin.cpp.function.CCFEMFunctionSpace

Bases: object

This class represents a cut and composite finite element function space (CCFEM) defined on one or more possibly intersecting meshes.

A CCFEM function space may be created from a set of standard function spaces by repeatedly calling add(), followed by a call to build(). Note that a CCFEM function space is not useful and its data structures are empty until build() has been called.

Create empty CCFEM function space

add()

Overloaded versions

  • add(function_space)

    Add function space (shared pointer version)

    Arguments
    function_space (FunctionSpace)

    The function space.

  • add(function_space)

    Add function space (reference version)

    Arguments
    function_space (FunctionSpace)

    The function space.

build()

Build CCFEM function space

collision_map_cut_cells()

Return the collision map for cut cells of the given part

Arguments
part (int)
The part number
Returns
std::map<unsigned int, std::vector<std::pair<std::size_t, unsigned int> > >
A map from cell indices of cut cells to a list of cutting cells. Each cutting cell is represented as a pair (part_number, cutting_cell_index).
covered_cells()

Return the list of covered cells for given part. The covered cells are defined as all cells that collide with the domain of any part with higher part number, but not with the boundary of that part; in other words cells that are completely covered by any other part (and which therefore are inactive).

Arguments
part (int)
The part number
Returns
numpy.array(int)
List of covered cell indices for given part
cut_cells()

Return the list of cut cells for given part. The cut cells are defined as all cells that collide with the boundary of any part with higher part number.

FIXME: Figure out whether this makes sense; a cell may collide with the boundary of part j but may still be covered completely by the domain of part j + 1. Possible solution is to for each part i check overlapping parts starting from the top and working back down to i + 1.

Arguments
part (int)
The part number
Returns
numpy.array(int)
List of cut cell indices for given part
dim()

Return dimension of the CCFEM function space

Returns
int
The dimension of the CCFEM function space.
dofmap()

Return CCFEM dofmap

Returns
CCFEMDofMap
The dofmap.
num_parts()

Return the number function spaces (parts) of the CCFEM function space

Returns
int
The number of function spaces (parts) of the CCFEM function space.
part()

Return function space (part) number i

Arguments
i (int)
The part number
Returns
FunctionSpace
Function space (part) number i
quadrature_rule_cut_cells()

Return quadrature rules for cut cells of the given part

Arguments
part (int)
The part number
Returns
std::map<unsigned int, std::pair<std::vector<double>, std::vector<double> > >
A map from cell indices of cut cells to a quadrature rules. Each quadrature rule is represented as a pair of an array of quadrature weights and a corresponding flattened array of quadrature points.
thisown

The membership flag

uncut_cells()

Return the list of uncut cells for given part. The uncut cells are defined as all cells that don’t collide with any cells in any other part with higher part number.

Arguments
part (int)
The part number
Returns
numpy.array(int)
List of uncut cell indices for given part