:py:mod:`symfem.finite_element` =============================== .. py:module:: symfem.finite_element .. autoapi-nested-parse:: Abstract finite element classes and functions. Module Contents --------------- Classes ~~~~~~~ .. autoapisummary:: symfem.finite_element.FiniteElement symfem.finite_element.CiarletElement symfem.finite_element.DirectElement symfem.finite_element.EnrichedElement symfem.finite_element.ElementBasisFunction Attributes ~~~~~~~~~~ .. autoapisummary:: symfem.finite_element.TabulatedBasis .. py:data:: TabulatedBasis .. py:exception:: NoTensorProduct Bases: :py:obj:`Exception` Error for element without a tensor representation. .. py:class:: FiniteElement(reference: symfem.references.Reference, order: int, space_dim: int, domain_dim: int, range_dim: int, range_shape: Optional[Tuple[int, Ellipsis]] = None) Bases: :py:obj:`abc.ABC` Abstract finite element. .. py:property:: maximum_degree :type: int :abstractmethod: Get the maximum degree of this polynomial set for the element. .. py:property:: name :type: str Get the name of the element. :returns: The name of the element's family .. py:attribute:: _float_basis_functions :type: Union[None, List[symfem.functions.AnyFunction]] .. py:attribute:: _value_scale :type: Union[None, sympy.core.expr.Expr] .. py:attribute:: names :type: List[str] :value: [] .. py:attribute:: references :type: List[str] :value: [] .. py:attribute:: last_updated .. py:attribute:: cache :value: True .. py:attribute:: _max_continuity_test_order :value: 4 .. py:method:: dof_plot_positions() -> List[symfem.geometry.PointType] :abstractmethod: Get the points to plot each DOF at on a DOF diagram. :returns: The DOF positions .. py:method:: dof_directions() -> List[Union[symfem.geometry.PointType, None]] :abstractmethod: Get the direction associated with each DOF. :returns: The DOF directions .. py:method:: dof_entities() -> List[Tuple[int, int]] :abstractmethod: Get the entities that each DOF is associated with. :returns: The entities .. py:method:: plot_dof_diagram(filename: Union[str, List[str]], plot_options: Dict[str, Any] = {}, **kwargs: Any) Plot a diagram showing the DOFs of the element. :param filename: The file name :param plot_options: Options for the plot :param kwargs: Keyword arguments .. py:method:: entity_dofs(entity_dim: int, entity_number: int) -> List[int] :abstractmethod: Get the numbers of the DOFs associated with the given entity. :param entity_dim: The dimension of the entity :param entity_number: The number of the entity :returns: The numbers of the DOFs associated with the entity .. py:method:: get_basis_functions(use_tensor_factorisation: bool = False) -> List[symfem.functions.AnyFunction] :abstractmethod: Get the basis functions of the element. :param use_tensor_factorisation: Should a tensor factorisation be used? :returns: The basis functions .. py:method:: get_basis_function(n: int) -> symfem.basis_functions.BasisFunction Get a single basis function of the element. :param n: The number of the basis function :returns: The basis function .. py:method:: tabulate_basis(points_in: symfem.geometry.SetOfPointsInput, order: str = 'xyzxyz') -> TabulatedBasis Evaluate the basis functions of the element at the given points. :param points_in: The points :param order: The order to return the values :returns: The tabulated basis functions .. py:method:: tabulate_basis_float(points_in: symfem.geometry.SetOfPointsInput) -> TabulatedBasis Evaluate the basis functions of the element at the given points in xyz,xyz order. :param points_in: The points :returns: The tabulated basis functions .. py:method:: plot_basis_function(n: int, filename: Union[str, List[str]], cell: Optional[symfem.references.Reference] = None, **kwargs: Any) Plot a diagram showing a basis function. :param n: The basis function number :param filename: The file name :param cell: The cell to push the basis function to and plot on :param kwargs: Keyword arguments .. py:method:: map_to_cell(vertices_in: symfem.geometry.SetOfPointsInput, basis: Optional[List[symfem.functions.AnyFunction]] = None, forward_map: Optional[symfem.geometry.PointType] = None, inverse_map: Optional[symfem.geometry.PointType] = None) -> List[symfem.functions.AnyFunction] :abstractmethod: Map the basis onto a cell using the appropriate mapping for the element. :param vertices_in: The vertices of the cell :param basis: The basis functions :param forward_map: The map from the reference to the cell :param inverse_map: The map to the reference from the cell :returns: The basis functions mapped to the cell .. py:method:: get_polynomial_basis() -> List[symfem.functions.AnyFunction] :abstractmethod: Get the symbolic polynomial basis for the element. :returns: The polynomial basis .. py:method:: test() Run tests for this element. .. py:method:: test_continuity() Test that this element has the correct continuity. .. py:method:: get_tensor_factorisation() -> List[Tuple[str, List[FiniteElement], List[int]]] Get the representation of the element as a tensor product. :returns: The tensor factorisation .. py:method:: _get_basis_functions_tensor() -> List[symfem.functions.AnyFunction] Compute the basis functions using the space's tensor product factorisation. :returns: The basis functions .. py:method:: init_kwargs() -> Dict[str, Any] Return the keyword arguments used to create this element. :returns: Keyword arguments dictionary .. py:class:: CiarletElement(reference: symfem.references.Reference, order: int, basis: List[symfem.functions.FunctionInput], dofs: symfem.functionals.ListOfFunctionals, domain_dim: int, range_dim: int, range_shape: Optional[Tuple[int, Ellipsis]] = None) Bases: :py:obj:`FiniteElement` Finite element defined using the Ciarlet definition. .. py:property:: maximum_degree :type: int Get the maximum degree of this polynomial set for the element. .. py:method:: entity_dofs(entity_dim: int, entity_number: int) -> List[int] Get the numbers of the DOFs associated with the given entity. :param entity_dim: The dimension of the entity :param entity_number: The number of the entity :returns: The numbers of the DOFs associated with the entity .. py:method:: dof_plot_positions() -> List[symfem.geometry.PointType] Get the points to plot each DOF at on a DOF diagram. :returns: The DOF positions .. py:method:: dof_directions() -> List[Union[symfem.geometry.PointType, None]] Get the direction associated with each DOF. :returns: The DOF directions .. py:method:: dof_entities() -> List[Tuple[int, int]] Get the entities that each DOF is associated with. :returns: The entities .. py:method:: get_polynomial_basis() -> List[symfem.functions.AnyFunction] Get the symbolic polynomial basis for the element. :returns: The polynomial basis .. py:method:: get_dual_matrix(inverse=False, caching=True) -> sympy.matrices.dense.MutableDenseMatrix Get the dual matrix. :param inverse: Should the dual matrix be inverted? :param caching: Should the result be cached :returns: The dual matrix .. py:method:: get_basis_functions(use_tensor_factorisation: bool = False) -> List[symfem.functions.AnyFunction] Get the basis functions of the element. :param use_tensor_factorisation: Should a tensor factorisation be used? :returns: The basis functions .. py:method:: plot_basis_function(n: int, filename: Union[str, List[str]], cell: Optional[symfem.references.Reference] = None, **kwargs: Any) Plot a diagram showing a basis function. :param n: The basis function number :param filename: The file name :param cell: The cell to push the basis function to and plot on :param kwargs: Keyword arguments .. py:method:: map_to_cell(vertices_in: symfem.geometry.SetOfPointsInput, basis: Optional[List[symfem.functions.AnyFunction]] = None, forward_map: Optional[symfem.geometry.PointType] = None, inverse_map: Optional[symfem.geometry.PointType] = None) -> List[symfem.functions.AnyFunction] Map the basis onto a cell using the appropriate mapping for the element. :param vertices_in: The vertices of the cell :param basis: The basis functions :param forward_map: The map from the reference to the cell :param inverse_map: The map to the reference from the cell :returns: The basis functions mapped to the cell .. py:method:: test() Run tests for this element. .. py:method:: test_dof_points() Test that DOF points are valid. .. py:method:: test_functional_entities() Test that the dof entities are valid and match the references of integrals. .. py:method:: test_functionals() Test that the functionals are satisfied by the basis functions. .. py:class:: DirectElement(reference: symfem.references.Reference, order: int, basis_functions: List[symfem.functions.FunctionInput], basis_entities: List[Tuple[int, int]], domain_dim: int, range_dim: int, range_shape: Optional[Tuple[int, Ellipsis]] = None) Bases: :py:obj:`FiniteElement` Finite element defined directly. .. py:property:: maximum_degree :type: int :abstractmethod: Get the maximum degree of this polynomial set for the element. .. py:attribute:: _basis_functions :type: List[symfem.functions.AnyFunction] .. py:method:: entity_dofs(entity_dim: int, entity_number: int) -> List[int] Get the numbers of the DOFs associated with the given entity. :param entity_dim: The dimension of the entity :param entity_number: The number of the entity :returns: The numbers of the DOFs associated with the entity .. py:method:: dof_plot_positions() -> List[symfem.geometry.PointType] Get the points to plot each DOF at on a DOF diagram. :returns: The DOF positions .. py:method:: dof_directions() -> List[Union[symfem.geometry.PointType, None]] Get the direction associated with each DOF. :returns: The DOF directions .. py:method:: dof_entities() -> List[Tuple[int, int]] Get the entities that each DOF is associated with. :returns: The entities .. py:method:: get_basis_functions(use_tensor_factorisation: bool = False) -> List[symfem.functions.AnyFunction] Get the basis functions of the element. :param use_tensor_factorisation: Should a tensor factorisation be used? :returns: The basis functions .. py:method:: map_to_cell(vertices_in: symfem.geometry.SetOfPointsInput, basis: Optional[List[symfem.functions.AnyFunction]] = None, forward_map: Optional[symfem.geometry.PointType] = None, inverse_map: Optional[symfem.geometry.PointType] = None) -> List[symfem.functions.AnyFunction] Map the basis onto a cell using the appropriate mapping for the element. :param vertices_in: The vertices of the cell :param basis: The basis functions :param forward_map: The map from the reference to the cell :param inverse_map: The map to the reference from the cell :returns: The basis functions mapped to the cell .. py:method:: get_polynomial_basis() -> List[symfem.functions.AnyFunction] :abstractmethod: Get the symbolic polynomial basis for the element. :returns: The polynomial basis .. py:method:: test() Run tests for this element. .. py:method:: test_independence() Test that the basis functions of this element are linearly independent. .. py:class:: EnrichedElement(subelements: List[FiniteElement]) Bases: :py:obj:`FiniteElement` Finite element defined directly. .. py:property:: maximum_degree :type: int :abstractmethod: Get the maximum degree of this polynomial set for the element. .. py:attribute:: _basis_functions :type: Optional[List[symfem.functions.AnyFunction]] .. py:method:: entity_dofs(entity_dim: int, entity_number: int) -> List[int] Get the numbers of the DOFs associated with the given entity. :param entity_dim: The dimension of the entity :param entity_number: The number of the entity :returns: The numbers of the DOFs associated with the entity .. py:method:: dof_plot_positions() -> List[symfem.geometry.PointType] Get the points to plot each DOF at on a DOF diagram. :returns: The DOF positions .. py:method:: dof_directions() -> List[Union[symfem.geometry.PointType, None]] Get the direction associated with each DOF. :returns: The DOF directions .. py:method:: dof_entities() -> List[Tuple[int, int]] Get the entities that each DOF is associated with. :returns: The entities .. py:method:: get_basis_functions(use_tensor_factorisation: bool = False) -> List[symfem.functions.AnyFunction] Get the basis functions of the element. :param use_tensor_factorisation: Should a tensor factorisation be used? :returns: The basis functions .. py:method:: map_to_cell(vertices_in: symfem.geometry.SetOfPointsInput, basis: Optional[List[symfem.functions.AnyFunction]] = None, forward_map: Optional[symfem.geometry.PointType] = None, inverse_map: Optional[symfem.geometry.PointType] = None) -> List[symfem.functions.AnyFunction] Map the basis onto a cell using the appropriate mapping for the element. :param vertices_in: The vertices of the cell :param basis: The basis functions :param forward_map: The map from the reference to the cell :param inverse_map: The map to the reference from the cell :returns: The basis functions mapped to the cell .. py:method:: get_polynomial_basis() -> List[symfem.functions.AnyFunction] :abstractmethod: Get the symbolic polynomial basis for the element. :returns: The polynomial basis .. py:method:: test() Run tests for this element. .. py:class:: ElementBasisFunction(element: FiniteElement, n: int) Bases: :py:obj:`symfem.basis_functions.BasisFunction` A basis function of a finite element. .. py:method:: get_function() -> symfem.functions.AnyFunction Get the actual basis function. :returns: The basis function