Detailed API description¶
This section contains API documentation for the most commonly used interfaces of the library.
Module: skfem.mesh¶
This module defines finite element meshes. Meshes can be created using built-in constructors or loaded from external formats using meshio. The supported types are
MeshTri
, triangular meshMeshQuad
, quadrilateral meshMeshTet
, tetrahedral meshMeshHex
, hexahedral meshMeshLine
, one-dimensional mesh
Default constructor creates a mesh for the unit square:
>>> from skfem.mesh import MeshTri
>>> MeshTri()
<skfem MeshTri1 object>
Number of elements: 2
Number of vertices: 4
Number of nodes: 4
Named boundaries [# facets]: left [1], bottom [1], right [1], top [1]
Each mesh type has several constructors; see the docstring, e.g.,
help(MeshTri)
or click MeshTri
in the online
documentation. Importing from external formats can be done with the
constructor load()
.
Abstract class: Mesh¶
- class skfem.mesh.Mesh(doflocs: ~numpy.ndarray, t: ~numpy.ndarray, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element.Element'>, affine: bool = False, sort_t: bool = False, validate: bool = True)¶
A mesh defined by
Element
class.Mesh
is defined as a combination of elements/cells by specifying the spatial locations of the finite element nodes.- elements_satisfying(test: Callable[[ndarray], ndarray]) ndarray ¶
Return elements whose midpoints satisfy some condition.
- Parameters:
test – A function which returns
True
for the element midpoints that are to be included in the return set.
- facets_satisfying(test: Callable[[ndarray], ndarray], boundaries_only: bool = False, normal: ndarray | None = None) ndarray ¶
Return facets whose midpoints satisfy some condition.
- Parameters:
test – A function which returns
True
for the facet midpoints that are to be included in the return set.boundaries_only – If
True
, include only boundary facets.normal – If given, used to orient the set of facets.
- classmethod load(filename: str, out: List[str] | None = None, **kwargs)¶
Load a mesh using meshio.
- Parameters:
filename – The filename of the mesh file.
out – Optional list of
meshio.Mesh
attribute names, overwrite with the corresponding data. E.g.,['point_data', 'cell_data']
.
- nodes_satisfying(test: Callable[[ndarray], ndarray], boundaries_only: bool = False) ndarray ¶
Return nodes that satisfy some condition.
- Parameters:
test – A function which returns
True
for the set of nodes that are to be included in the return set.boundaries_only – If
True
, include only boundary facets.
- refined(times_or_ix: int | ndarray = 1)¶
Return a refined mesh.
- Parameters:
times_or_ix – Either an integer giving the number of uniform refinements or an array of element indices for adaptive refinement.
- save(filename: str | PathLike, point_data: Dict[str, ndarray] | None = None, cell_data: Dict[str, ndarray] | None = None, **kwargs) None ¶
Export the mesh and fields using meshio.
- Parameters:
filename – The output filename, with suffix determining format; e.g. .msh, .vtk, .xdmf
point_data – Data related to the vertices of the mesh.
cell_data – Data related to the elements of the mesh.
Class: MeshTri¶
- class skfem.mesh.MeshTri1(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_tri.element_tri_p1.ElementTriP1'>, affine: bool = True, sort_t: bool = True, validate: bool = True)¶
A standard first-order triangular mesh.
- __init__(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_tri.element_tri_p1.ElementTriP1'>, affine: bool = True, sort_t: bool = True, validate: bool = True) None ¶
- classmethod init_circle(nrefs: int = 3, smoothed: bool = False) Mesh2D ¶
Initialize a circle mesh.
Works by repeatedly refining the following mesh and moving new nodes to the boundary:
* / | \ / | \ / | \ *------O------* \ | / \ | / \ | / *
- Parameters:
nrefs – Number of refinements, by default 3.
smoothed – If
True
, apply smoothing after each refine.
- classmethod init_lshaped() Mesh2D ¶
Initialize a mesh for the L-shaped domain.
The mesh topology is as follows:
*-------* | \ | | \ | | \ | |-------O-------* | / | \ | | / | \ | | / | \ | *---------------*
- classmethod init_refdom()¶
Initialize a mesh corresponding to the reference domain.
- classmethod init_sqsymmetric() Mesh2D ¶
Initialize a symmetric mesh of the unit square.
The mesh topology is as follows:
*------*------* |\ | /| | \ | / | | \ | / | *------*------* | / | \ | | / | \ | |/ | \| O------*------*
- classmethod init_symmetric() Mesh2D ¶
Initialize a symmetric mesh of the unit square.
The mesh topology is as follows:
*------------* |\ /| | \ / | | \ / | | * | | / \ | | / \ | |/ \| O------------*
- classmethod init_tensor(x: ndarray, y: ndarray)¶
Initialize a tensor product mesh.
The mesh topology is as follows:
*---------------* |'-.|'-.|`'---._| |---+---+-------| |\ |\ |'. | | \ | \ | '-. | | \| \| '.| *---------------*
- Parameters:
x – The nodal coordinates in dimension x.
y – The nodal coordinates in dimension y.
- classmethod load(filename: str, out: List[str] | None = None, **kwargs)¶
Load a mesh using meshio.
- Parameters:
filename – The filename of the mesh file.
out – Optional list of
meshio.Mesh
attribute names, overwrite with the corresponding data. E.g.,['point_data', 'cell_data']
.
- class skfem.mesh.MeshTri2(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_tri.element_tri_p2.ElementTriP2'>, affine: bool = False, sort_t: bool = False, validate: bool = True)¶
A quadratic triangular mesh.
- classmethod init_circle(nrefs: int = 3, smoothed: bool = False) MeshTri2 ¶
Initialize a circle mesh.
Works by repeatedly refining the following mesh and moving new nodes to the boundary:
* / | \ / | \ / | \ *------O------* \ | / \ | / \ | / *
- Parameters:
nrefs – Number of refinements, by default 3.
smoothed – If
True
, apply smoothing after each refine.
- classmethod load(filename: str, out: List[str] | None = None, **kwargs)¶
Load a mesh using meshio.
- Parameters:
filename – The filename of the mesh file.
out – Optional list of
meshio.Mesh
attribute names, overwrite with the corresponding data. E.g.,['point_data', 'cell_data']
.
Class: MeshQuad¶
- class skfem.mesh.MeshQuad1(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_quad.element_quad1.ElementQuad1'>, affine: bool = False, sort_t: bool = False, validate: bool = True)¶
A standard first-order quadrilateral mesh.
If
t
is provided, order of vertices in each element should match the numbering:3---2 | | 0---1
- __init__(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_quad.element_quad1.ElementQuad1'>, affine: bool = False, sort_t: bool = False, validate: bool = True) None ¶
- classmethod init_refdom()¶
Initialize a mesh corresponding to the reference domain.
- classmethod init_tensor(x: ndarray, y: ndarray)¶
Initialize a tensor product mesh.
The mesh topology is as follows:
*-------------* | | | | |---+--+------| | | | | | | | | | | | | *-------------*
- Parameters:
x – The nodal coordinates in dimension x.
y – The nodal coordinates in dimension y.
- classmethod load(filename: str, out: List[str] | None = None, **kwargs)¶
Load a mesh using meshio.
- Parameters:
filename – The filename of the mesh file.
out – Optional list of
meshio.Mesh
attribute names, overwrite with the corresponding data. E.g.,['point_data', 'cell_data']
.
- to_meshtri(x: ndarray | None = None, style: str | None = None)¶
Split each quadrilateral into two triangles.
- Parameters:
x – Elementwise constant function to preserve. If given, returned as the second additional output parameter.
style – Optionally, specify a splitting style ‘x’ for crisscross splitting.
- class skfem.mesh.MeshQuad2(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_quad.element_quad2.ElementQuad2'>, affine: bool = False, sort_t: bool = False, validate: bool = True)¶
A quadratic quadrilateral mesh.
- classmethod load(filename: str, out: List[str] | None = None, **kwargs)¶
Load a mesh using meshio.
- Parameters:
filename – The filename of the mesh file.
out – Optional list of
meshio.Mesh
attribute names, overwrite with the corresponding data. E.g.,['point_data', 'cell_data']
.
Class: MeshTet¶
- class skfem.mesh.MeshTet1(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_tet.element_tet_p1.ElementTetP1'>, affine: bool = True, sort_t: bool = False, validate: bool = True)¶
A standard first-order tetrahedral mesh.
- __init__(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_tet.element_tet_p1.ElementTetP1'>, affine: bool = True, sort_t: bool = False, validate: bool = True) None ¶
- classmethod init_ball(nrefs: int = 3)¶
Initialize a ball mesh.
- Parameters:
nrefs – Number of refinements, by default 3.
- classmethod init_refdom()¶
Initialize a mesh corresponding to the reference domain.
- classmethod init_tensor(x: ndarray, y: ndarray, z: ndarray)¶
Initialize a tensor product mesh.
- Parameters:
x – The nodal coordinates in dimension x.
y – The nodal coordinates in dimension y.
z – The nodal coordinates in dimension z.
- classmethod load(filename: str, out: List[str] | None = None, **kwargs)¶
Load a mesh using meshio.
- Parameters:
filename – The filename of the mesh file.
out – Optional list of
meshio.Mesh
attribute names, overwrite with the corresponding data. E.g.,['point_data', 'cell_data']
.
- class skfem.mesh.MeshTet2(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_tet.element_tet_p2.ElementTetP2'>, affine: bool = False, sort_t: bool = False, validate: bool = True)¶
A quadratic tetrahedral mesh.
- classmethod init_ball(nrefs: int = 3) MeshTet2 ¶
Initialize a ball mesh.
- Parameters:
nrefs – Number of refinements, by default 3.
- classmethod load(filename: str, out: List[str] | None = None, **kwargs)¶
Load a mesh using meshio.
- Parameters:
filename – The filename of the mesh file.
out – Optional list of
meshio.Mesh
attribute names, overwrite with the corresponding data. E.g.,['point_data', 'cell_data']
.
Class: MeshHex¶
- class skfem.mesh.MeshHex1(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_hex.element_hex1.ElementHex1'>, affine: bool = False, sort_t: bool = False, validate: bool = True)¶
A standard first-order hexahedral mesh.
If
t
is provided, order of vertices in each element should match the numbering:2---6 / /| 4---7 3 | |/ 1---5
- __init__(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_hex.element_hex1.ElementHex1'>, affine: bool = False, sort_t: bool = False, validate: bool = True) None ¶
- classmethod init_tensor(x: ndarray, y: ndarray, z: ndarray)¶
Initialize a tensor product mesh.
- Parameters:
x – The nodal coordinates in dimension x.
y – The nodal coordinates in dimension y.
z – The nodal coordinates in dimension z.
- classmethod load(filename: str, out: List[str] | None = None, **kwargs)¶
Load a mesh using meshio.
- Parameters:
filename – The filename of the mesh file.
out – Optional list of
meshio.Mesh
attribute names, overwrite with the corresponding data. E.g.,['point_data', 'cell_data']
.
- to_meshtet()¶
Split each hexahedron into six tetrahedra.
Class: MeshLine¶
- class skfem.mesh.MeshLine1(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_line.element_line_p1.ElementLineP1'>, affine: bool = True, sort_t: bool = False, validate: bool = True)¶
A one-dimensional mesh.
- __init__(doflocs: ~numpy.ndarray = <factory>, t: ~numpy.ndarray = <factory>, _boundaries: ~typing.Dict[str, ~numpy.ndarray] | None = None, _subdomains: ~typing.Dict[str, ~numpy.ndarray] | None = None, elem: ~typing.Type[~skfem.element.element.Element] = <class 'skfem.element.element_line.element_line_p1.ElementLineP1'>, affine: bool = True, sort_t: bool = False, validate: bool = True) None ¶
Module: skfem.assembly¶
This module performs the finite element assembly. The basic workflow is the following:
Initialize
Mesh
andElement
.
>>> import skfem as fem
>>> m = fem.MeshTri()
>>> e = fem.ElementTriP1()
>>> m
<skfem MeshTri1 object>
Number of elements: 2
Number of vertices: 4
Number of nodes: 4
Named boundaries [# facets]: left [1], bottom [1], right [1], top [1]
Create
CellBasis
orFacetBasis
objects.
>>> basis = fem.CellBasis(m, e)
Define the forms using
BilinearForm
,LinearForm
, orFunctional
.
>>> form_a = fem.BilinearForm(lambda u, v, w: u * v)
>>> form_l = fem.LinearForm(lambda v, w: w.x[0] ** 2 * v)
Mathematically the above forms are
Create the matrices/vectors using
assemble()
.
>>> A = form_a.assemble(basis)
>>> b = form_l.assemble(basis)
>>> A.toarray()
array([[0.08333333, 0.04166667, 0.04166667, 0. ],
[0.04166667, 0.16666667, 0.08333333, 0.04166667],
[0.04166667, 0.08333333, 0.16666667, 0.04166667],
[0. , 0.04166667, 0.04166667, 0.08333333]])
>>> b
array([0.0162037 , 0.15046296, 0.06712963, 0.09953704])
- skfem.assembly.asm(form: ~skfem.assembly.form.form.Form, *args, to=<function _sum>, **kwargs) Any ¶
Perform finite element assembly.
A shorthand for
skfem.assembly.Form.assemble()
which, in addition, supports assembling multiple bases at once and summing the result.
Abstract class: AbstractBasis¶
Subclasses of AbstractBasis
represent a global
finite element basis evaluated at quadrature points.
- class skfem.assembly.basis.AbstractBasis(mesh: ~skfem.mesh.mesh.Mesh, elem: ~skfem.element.element.Element, mapping: ~skfem.mapping.mapping.Mapping | None = None, intorder: int | None = None, quadrature: ~typing.Tuple[~numpy.ndarray, ~numpy.ndarray] | None = None, refdom: ~typing.Type[~skfem.refdom.Refdom] = <class 'skfem.refdom.Refdom'>, dofs: ~skfem.assembly.dofs.Dofs | None = None, disable_doflocs: bool = False)¶
Finite element basis at global quadrature points.
Please see the following implementations:
CellBasis
, basis functions inside elementsFacetBasis
, basis functions on boundaryInteriorFacetBasis
, basis functions on facets inside the domain
- get_dofs(facets: Any | None = None, elements: Any | None = None, nodes: Any | None = None, skip: List[str] | None = None) Any ¶
Find global DOF numbers.
Accepts an array of facet/element indices. However, various argument types can be turned into an array of facet/element indices.
Get all boundary DOFs:
>>> import numpy as np >>> from skfem import MeshTri, Basis, ElementTriP1 >>> m = MeshTri().refined() >>> basis = Basis(m, ElementTriP1()) >>> basis.get_dofs().flatten() array([0, 1, 2, 3, 4, 5, 7, 8])
Get DOFs via a function query:
>>> import numpy as np >>> from skfem import MeshTri, Basis, ElementTriP1 >>> m = MeshTri().refined() >>> basis = Basis(m, ElementTriP1()) >>> basis.get_dofs(lambda x: np.isclose(x[0], 0)).flatten() array([0, 2, 5])
Get DOFs via named boundaries:
>>> import numpy as np >>> from skfem import MeshTri, Basis, ElementTriP1 >>> m = (MeshTri() ... .refined() ... .with_boundaries({'left': lambda x: np.isclose(x[0], 0)})) >>> basis = Basis(m, ElementTriP1()) >>> basis.get_dofs('left').flatten() array([0, 2, 5])
Get DOFs via named subdomains:
>>> from skfem import MeshTri, Basis, ElementTriP1 >>> m = (MeshTri() ... .refined() ... .with_subdomains({'left': lambda x: x[0] < .5})) >>> basis = Basis(m, ElementTriP1()) >>> basis.get_dofs(elements='left').flatten() array([0, 2, 4, 5, 6, 8])
Further reduce the set of DOFs:
>>> import numpy as np >>> from skfem import MeshTri, Basis, ElementTriArgyris >>> m = MeshTri().refined() >>> basis = Basis(m, ElementTriArgyris()) >>> basis.get_dofs(lambda x: np.isclose(x[0], 0)).nodal.keys() dict_keys(['u', 'u_x', 'u_y', 'u_xx', 'u_xy', 'u_yy']) >>> basis.get_dofs(lambda x: np.isclose(x[0], 0)).all(['u', 'u_x']) array([ 0, 1, 12, 13, 30, 31])
Skip some DOF names altogether:
>>> import numpy as np >>> from skfem import MeshTri, Basis, ElementTriArgyris >>> m = MeshTri().refined() >>> basis = Basis(m, ElementTriArgyris()) >>> basis.get_dofs(lambda x: np.isclose(x[0], 0), ... skip=['u_x', 'u_y']).nodal.keys() dict_keys(['u', 'u_xx', 'u_xy', 'u_yy'])
Combine several boundaries into one call:
>>> import numpy as np >>> from skfem import MeshTri, Basis, ElementTriP1 >>> m = (MeshTri() ... .with_boundaries({'left': lambda x: np.isclose(x[0], 0), ... 'right': lambda x: np.isclose(x[0], 1)})) >>> basis = Basis(m, ElementTriP1()) >>> basis.get_dofs({'left', 'right'}).flatten() array([0, 1, 2, 3])
- Parameters:
facets – An array of facet indices. If
None
, find facets byself.mesh.boundary_facets()
. If callable, callself.mesh.facets_satisfying(facets)
to get the facets. If string, useself.mesh.boundaries[facets]
to get the facets. If list, tuple or set, use the combined facet indices.elements – An array of element indices. See above.
nodes – An array of node indices.
skip – List of dofnames to skip.
Class: CellBasis¶
- class skfem.assembly.CellBasis(mesh: Mesh, elem: Element, mapping: Mapping | None = None, intorder: int | None = None, elements: Any | None = None, quadrature: Tuple[ndarray, ndarray] | None = None, dofs: Dofs | None = None, disable_doflocs: bool = False)¶
For fields defined inside the domain.
CellBasis
object is a combination ofMesh
andElement
.>>> from skfem import * >>> m = MeshTri.init_symmetric() >>> e = ElementTriP1() >>> basis = CellBasis(m, e)
The resulting objects are used in the assembly.
>>> from skfem.models.poisson import laplace >>> K = asm(laplace, basis) >>> K.shape (5, 5)
- __init__(mesh: Mesh, elem: Element, mapping: Mapping | None = None, intorder: int | None = None, elements: Any | None = None, quadrature: Tuple[ndarray, ndarray] | None = None, dofs: Dofs | None = None, disable_doflocs: bool = False)¶
Combine
Mesh
andElement
into a set of precomputed global basis functions.- Parameters:
mesh – An object of type
Mesh
.elem – An object of type
Element
.mapping – An object of type
skfem.mapping.Mapping
. If None, uses mesh.mapping.intorder – Optional integration order, i.e. the degree of polynomials that are integrated exactly by the used quadrature. Not used if quadrature is specified.
elements – Optional subset of element indices.
quadrature – Optional tuple of quadrature points and weights.
dofs – Optional
Dofs
object.disable_doflocs – If True, the computation of global DOF locations is disabled. This may save memory on large meshes if DOF locations are not required.
- interpolate(w: ndarray) DiscreteField | Tuple[DiscreteField, ...] ¶
Interpolate a solution vector to quadrature points.
Useful when a solution vector is needed in the forms, e.g., when evaluating functionals or when solving nonlinear problems.
- Parameters:
w – A solution vector.
- project(interp, elements=None, dtype=None)¶
Perform \(L^2\) projection onto the basis.
See Performing projections for more information.
- Parameters:
interp – An object of type
DiscreteField
which is a function (to be projected) evaluated at global quadrature points. If a function is given, thenDiscreteField
is created by passing an array of global quadrature point locations to the function.elements – Optionally perform the projection on a subset of elements. The values of the remaining DOFs are zero.
dtype – Set to np.complex64 or similar to use complex numbers.
Class: FacetBasis¶
- skfem.assembly.BoundaryFacetBasis¶
alias of
FacetBasis
- class skfem.assembly.FacetBasis(mesh: Mesh, elem: Element, mapping: Mapping | None = None, intorder: int | None = None, quadrature: Tuple[ndarray, ndarray] | None = None, facets: Any | None = None, dofs: Dofs | None = None, side: int = 0, disable_doflocs: bool = False)¶
For integrating over facets of the mesh. Usually over the boundary.
- __init__(mesh: Mesh, elem: Element, mapping: Mapping | None = None, intorder: int | None = None, quadrature: Tuple[ndarray, ndarray] | None = None, facets: Any | None = None, dofs: Dofs | None = None, side: int = 0, disable_doflocs: bool = False)¶
Precomputed global basis on boundary facets.
- Parameters:
mesh – An object of type
Mesh
.elem – An object of type
Element
.mapping – An object of type
skfem.mapping.Mapping
. If None, uses mesh.mapping.intorder – Optional integration order, i.e. the degree of polynomials that are integrated exactly by the used quadrature. Not used if quadrature is specified.
quadrature – Optional tuple of quadrature points and weights.
facets – Optional subset of facet indices.
dofs – Optional
Dofs
object.disable_doflocs – If True, the computation of global DOF locations is disabled. This may save memory on large meshes if DOF locations are not required.
Class: InteriorFacetBasis¶
- class skfem.assembly.InteriorFacetBasis(mesh: Mesh, elem: Element, mapping: Mapping | None = None, intorder: int | None = None, quadrature: Tuple[ndarray, ndarray] | None = None, facets: Any | None = None, dofs: Dofs | None = None, side: int = 0, disable_doflocs: bool = False)¶
For evaluating integrals over interior facets.
Useful for, e.g., a posteriori error estimators or implementing interior penalty/discontinuous Galerkin methods.
Abstract class: Form¶
Class: BilinearForm¶
- class skfem.assembly.BilinearForm(form: ~typing.Callable | ~skfem.assembly.form.form.Form | None = None, dtype: ~typing.Type[~numpy.float64] | ~typing.Type[~numpy.complex64] = <class 'numpy.float64'>, nthreads: int = 0, **params)¶
A bilinear form for finite element assembly.
Bilinear forms are defined using functions that takes three arguments: trial function
u
, test functionv
, and a dictionary of additional parametersw
.>>> from skfem import BilinearForm, Basis, MeshTri, ElementTriP1 >>> form = BilinearForm(lambda u, v, _: u * v) >>> form.assemble(Basis(MeshTri(), ElementTriP1())).toarray() array([[0.08333333, 0.04166667, 0.04166667, 0. ], [0.04166667, 0.16666667, 0.08333333, 0.04166667], [0.04166667, 0.08333333, 0.16666667, 0.04166667], [0. , 0.04166667, 0.04166667, 0.08333333]])
Alternatively, you can use
BilinearForm
as a decorator:>>> @BilinearForm ... def form(u, v, _): ... return u * v
Inside the form definition,
u
andv
are arrays containing the basis function values at quadrature points. They also contain the values of the derivatives as attributes:>>> @BilinearForm ... def form(u, v, _): ... # u.grad[0] is first derivative with respect to x, and so on ... return u.grad[0] * v.grad[0] + u.grad[1] * v.grad[1] # laplacian
We suggest you to use helper functions from
skfem.helpers
to make the forms more readable:>>> from skfem.helpers import dot, grad >>> @BilinearForm ... def form(u, v, w): ... return dot(grad(u), grad(v))
Class: LinearForm¶
- class skfem.assembly.LinearForm(form: ~typing.Callable | ~skfem.assembly.form.form.Form | None = None, dtype: ~typing.Type[~numpy.float64] | ~typing.Type[~numpy.complex64] = <class 'numpy.float64'>, nthreads: int = 0, **params)¶
A linear form for finite element assembly.
Used similarly as
BilinearForm
with the expection that forms take two parametersv
andw
.
Class: Functional¶
- class skfem.assembly.Functional(form: ~typing.Callable | ~skfem.assembly.form.form.Form | None = None, dtype: ~typing.Type[~numpy.float64] | ~typing.Type[~numpy.complex64] = <class 'numpy.float64'>, nthreads: int = 0, **params)¶
A functional for postprocessing finite element solution.
Used similarly as
BilinearForm
with the expection that forms take one parameterw
.- elemental(v: AbstractBasis, **kwargs) ndarray ¶
Evaluate the functional elementwise.
Module: skfem.element¶
The module skfem.element
defines finite elements in a very generic
sense.
In order to use an element, you simply initialize the respective object and
pass it to the constructor of CellBasis
or
FacetBasis
. See below for a list of supported
elements.
|
\(H^1\)-conforming basis defined through a reference element. |
|
Use the same element for each dimension. |
|
\(H(div)\)-conforming basis defined through a reference element. |
|
\(H(curl)\)-conforming basis defined through a reference element. |
|
Elements defined implicitly through global degrees-of-freedom. |
|
Turn a finite element discontinuous by cutting connectivity. |
|
Combine multiple elements. |
|
Piecewise linear element. |
|
Piecewise quadratic element. |
|
Piecewise cubic element. |
|
Piecewise quartic element. |
|
Piecewise constant element. |
|
The nonconforming Crouzeix-Raviart element. |
|
alias of |
|
alias of |
|
The Morley element. |
|
15-parameter nonconforming plate element. |
|
The Argyris element. |
|
alias of |
|
The Hermite element with 3 DOFs per vertex and one interior DOF. |
|
Constant element for the mesh skeleton. |
|
Linear element for the mesh skeleton. |
|
The lowest order Brezzi-Douglas-Marini element. |
|
Piecewise constant element. |
|
Bilinear element. |
|
Biquadratic element. |
|
Serendipity element with 8 DOFs. |
|
Piecewise \(p\)'th order element. |
|
The Bogner-Fox-Schmit element. |
|
alias of |
|
Piecewise constant element. |
|
Piecewise linear element. |
|
Piecewise quadratic element. |
|
alias of |
|
alias of |
|
The MINI element, i.e. linear element with additional bubble DOF. |
|
Nonconforming Crouzeix-Raviart element. |
|
Conforming Crouzeix-Raviart element. |
|
Piecewise constant element. |
|
Trilinear element. |
|
Triquadratic element. |
|
Serendipity element with 20 DOFs. |
|
Piecewise constant element. |
|
Piecewise linear element. |
|
Piecewise quadratic element. |
|
Piecewise \(p\)'th order element. |
|
\(H^2\)-conforming element with 4 DOFs. |
|
The MINI element, i.e. linear element with additional bubble DOF. |
Note
The element global basis is calculated at quadrature points and stored inside
DiscreteField
objects.
The different element types precalculate different fields of
DiscreteField
. E.g., for \(H(div)\)
finite elements it is natural to precalculate DiscreteField.div
.
The high order derivatives are created only when using subclasses of
ElementGlobal
.
- class skfem.element.DiscreteField(value=None, grad=None, div=None, curl=None, hess=None, grad3=None, grad4=None, grad5=None, grad6=None)¶
Module: skfem.utils¶
Function: solve¶
- skfem.utils.solve(A: spmatrix, b: ndarray | spmatrix, x: ndarray | None = None, I: ndarray | None = None, solver: Callable[[...], ndarray] | Callable[[...], Tuple[ndarray, ndarray]] | None = None, **kwargs) ndarray | Tuple[ndarray, ndarray] ¶
Solve a linear system or a generalized eigenvalue problem.
The remaining keyword arguments are passed to the solver.
- Parameters:
A – The system matrix
b – The right hand side vector or the mass matrix of a generalized eigenvalue problem.
solver – Choose one of the following solvers:
skfem.utils.solver_direct_scipy()
(default),skfem.utils.solver_eigen_scipy()
(default),skfem.utils.solver_iter_pcg()
,skfem.utils.solver_iter_krylov()
.
Function: condense¶
- skfem.utils.condense(A: spmatrix, b: ndarray | spmatrix | None = None, x: ndarray | None = None, I: ndarray | DofsView | Dict[str, DofsView] | None = None, D: ndarray | DofsView | Dict[str, DofsView] | None = None, expand: bool = True) spmatrix | Tuple[spmatrix, ndarray] | Tuple[spmatrix, spmatrix] | Tuple[spmatrix, ndarray, ndarray] | Tuple[spmatrix, ndarray, ndarray, ndarray] | Tuple[spmatrix, spmatrix, ndarray, ndarray] ¶
Eliminate degrees-of-freedom from a linear system.
The user should provide the linear system
A
andb
and either the set of DOFs to eliminate (D
) or the set of DOFs to keep (I
). Optionally, nonzero values for the eliminated DOFs can be supplied viax
.Note
Supports also generalized eigenvalue problems where
b
is a matrix.Example
Suppose that the solution vector \(x\) can be split as
\[\begin{split}x = \begin{bmatrix} x_I\\ x_D \end{bmatrix}\end{split}\]where \(x_D\) are known and \(x_I\) are unknown. This allows splitting the linear system as
\[\begin{split}\begin{bmatrix} A_{II} & A_{ID}\\ A_{DI} & A_{DD} \end{bmatrix} \begin{bmatrix} x_I\\ x_D \end{bmatrix} = \begin{bmatrix} b_I\\ b_D \end{bmatrix}\end{split}\]which leads to the condensed system
\[A_{II} x_I = b_I - A_{ID} x_D.\]As an example, let us assemble the matrix \(A\) and the vector \(b\) corresponding to the Poisson equation \(-\Delta u = 1\).
>>> import skfem as fem >>> from skfem.models.poisson import laplace, unit_load >>> m = fem.MeshTri().refined(2) >>> basis = fem.CellBasis(m, fem.ElementTriP1()) >>> A = laplace.assemble(basis) >>> b = unit_load.assemble(basis)
The condensed system is obtained with
skfem.utils.condense()
. Below we provide the DOFs to eliminate via the keyword argumentD
.>>> AII, bI, xI, I = fem.condense(A, b, D=m.boundary_nodes()) >>> AII.toarray() array([[ 4., 0., 0., 0., -1., -1., -1., -1., 0.], [ 0., 4., 0., 0., -1., 0., -1., 0., 0.], [ 0., 0., 4., 0., 0., -1., 0., -1., 0.], [ 0., 0., 0., 4., -1., -1., 0., 0., 0.], [-1., -1., 0., -1., 4., 0., 0., 0., 0.], [-1., 0., -1., -1., 0., 4., 0., 0., 0.], [-1., -1., 0., 0., 0., 0., 4., 0., -1.], [-1., 0., -1., 0., 0., 0., 0., 4., -1.], [ 0., 0., 0., 0., 0., 0., -1., -1., 4.]]) >>> bI array([0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625])
By default, the eliminated DOFs are set to zero. Different values can be provided through the keyword argument
x
; see Example 14: Laplace with inhomogeneous boundary conditions.- Parameters:
A – The system matrix
b – The right hand side vector, or zero if x is given, or the mass matrix for generalized eigenvalue problems.
x – The values of the condensed degrees-of-freedom. If not given, assumed to be zero.
I – The set of degree-of-freedom indices to include.
D – The set of degree-of-freedom indices to dismiss.
expand – If
True
(default), returns also x and I. As a consequence,skfem.utils.solve()
will expand the solution vector automatically.
- Returns:
The condensed linear system and (optionally) information about the boundary values.
- Return type:
CondensedSystem
Function: enforce¶
- skfem.utils.enforce(A: spmatrix, b: ndarray | spmatrix | None = None, x: ndarray | None = None, I: ndarray | DofsView | Dict[str, DofsView] | None = None, D: ndarray | DofsView | Dict[str, DofsView] | None = None, diag: float = 1.0, overwrite: bool = False) spmatrix | Tuple[spmatrix, ndarray] | Tuple[spmatrix, spmatrix] ¶
Enforce degrees-of-freedom of a linear system.
An alternative to
condense()
which sets the matrix diagonals to one and right-hand side vector to the enforced degree-of-freedom value.Note
The original system is both returned (for compatibility with
skfem.utils.solve()
) and optionally (if overwrite) modified (for performance).- Parameters:
A – The system matrix
b – Optionally, the right hand side vector.
x – The values of the enforced degrees-of-freedom. If not given, assumed to be zero.
I – Specify either this or
D
: The set of degree-of-freedom indices to solve for.D – Specify either this or
I
: The set of degree-of-freedom indices to enforce (rows/diagonal set to zero/one).overwrite – Optionally, the original system is both modified (for performance) and returned (for compatibility with
skfem.utils.solve()
). By default,False
.
- Returns:
A linear system with the enforced rows/diagonals set to zero/one.
- Return type:
LinearSystem
Module: skfem.helpers¶
Helper functions for defining forms.
- skfem.helpers.grad(u: DiscreteField)¶
Gradient.
- skfem.helpers.div(u: DiscreteField)¶
Divergence.
- skfem.helpers.curl(u: DiscreteField)¶
Curl.
- skfem.helpers.d(u: DiscreteField)¶
Gradient, divergence or curl.
- skfem.helpers.dd(u: DiscreteField)¶
Hessian (for
ElementGlobal
).
- skfem.helpers.ddd(u: DiscreteField)¶
Third derivative (for
ElementGlobal
).
- skfem.helpers.dddd(u: DiscreteField)¶
Fourth derivative (for
ElementGlobal
).
- skfem.helpers.sym_grad(u: DiscreteField)¶
Symmetric gradient.
- skfem.helpers.dot(u: DiscreteField | ndarray, v: DiscreteField | ndarray)¶
Dot product.
- skfem.helpers.ddot(u: DiscreteField | ndarray, v: DiscreteField | ndarray)¶
Double dot product.
- skfem.helpers.dddot(u: DiscreteField | ndarray, v: DiscreteField | ndarray)¶
Triple dot product.
- skfem.helpers.mul(A: DiscreteField | ndarray, x: DiscreteField | ndarray)¶
Matrix multiplication.
- skfem.helpers.trace(T)¶
Trace of matrix.
- skfem.helpers.transpose(T)¶
Transpose of matrix.
- skfem.helpers.prod(u: DiscreteField | ndarray, v: DiscreteField | ndarray, w: DiscreteField | ndarray | None = None)¶
Tensor product.
- skfem.helpers.inv(A)¶
Inverse of an array A over trailing axis (if any).