scikit-fem 2.2.0

Simple finite element assemblers

scikit-fem is a lightweight Python 3.6+ library for performing finite element assembly. Its main purpose is the transformation of bilinear forms into sparse matrices and linear forms into vectors. The library supports triangular, quadrilateral, tetrahedral and hexahedral meshes as well as one-dimensional problems.

The library fills an important gap in the spectrum of finite element codes. The library is lightweight meaning that it has minimal dependencies. It contains no compiled code meaning that it's easy to install and use on all platforms that support NumPy. Despite being fully interpreted, the code has a reasonably good performance.

The following benchmark (docs/examples/performance.py) demonstrates the time spent on finite element assembly in comparison to the time spent on linear solve. The given numbers were calculated using a ThinkPad X1 Carbon laptop (7th gen). Note that the timings are only illustrative as they depend on, e.g., the type of element used, the number of quadrature points used, the type of linear solver, and the complexity of the forms. This benchmark solves the Laplace equation using linear tetrahedral elements and the default direct sparse solver of scipy.sparse.linalg.spsolve.

Degrees-of-freedom	Assembly (s)	Linear solve (s)
64	0.00155	0.00073
125	0.00203	0.00072
216	0.00276	0.00081
512	0.00589	0.00127
1000	0.01076	0.00247
1728	0.02063	0.00538
4096	0.04805	0.04241
8000	0.09804	0.16269
15625	0.20347	0.87741
32768	0.46399	5.98163
64000	1.00143	36.47855
125000	2.05274	nan
262144	4.48825	nan
512000	8.82814	nan
1030301	18.25461	nan

Examples

Forms are defined using an intuitive syntax:

from skfem import BilinearForm
from skfem.helpers import dot, grad

@BilinearForm
def laplace(u, v, w):
    return dot(grad(u), grad(v))

Meshes can be initialized manually, loaded from external files using meshio, or created with the help of special constructors:

import numpy as np
from skfem import MeshLine, MeshTri, MeshTet

mesh = MeshLine(np.array([0.0, 0.5, 1.0]))
mesh = MeshTri.load("docs/examples/square.msh")
mesh = MeshTet.init_tensor(*((np.linspace(0, 1, 60),) * 3))

We support many common finite elements. Below the stiffness matrix is assembled using second-order tetrahedra:

from skfem import InteriorBasis, ElementTetP2

basis = InteriorBasis(mesh, ElementTetP2())
A = laplace.assemble(basis)  # type: scipy.sparse.csr_matrix

More examples can be found in the gallery.

Documentation

The project is documented using Sphinx. A recent version of the documentation can be found from Read the Docs.

Getting help

If you encounter an issue and cannot find help from the documentation, you can use the Github issue tracker to ask questions. Try to provide a snippet of code which fails and include also the version of the library you are using. The version can be found as follows:

python -c "import pkg_resources; print(pkg_resources.get_distribution('scikit-fem').version)"

Installation

The most recent release can be installed simply by pip install scikit-fem. For more cutting edge features, you can clone this repository.

Dependencies

The minimal dependencies for installing scikit-fem are numpy, scipy and meshio. In addition, many examples use matplotlib for visualization. Some examples demonstrate the use of other external packages; see our CI job definition for a full list of test dependencies.

Testing

The tests are run by Github Actions (see .github/workflows/main.yml).

The Makefile in the repository root has instructions for running the testing container locally using docker. For example, use make test_py38 to run the tests using py38 branch from kinnala/scikit-fem-docker-action.

Licensing

The contents of skfem/ and the PyPI package scikit-fem are licensed under the 3-clause BSD license. Some examples under docs/examples/ have a different license, see LICENSE.md for more information.

Acknowledgements

This project was started while working under a grant from the Finnish Cultural Foundation. The approach used in the finite element assembly has been inspired by the work of A. Hannukainen and M. Juntunen.

Contributing

We are happy to welcome any contributions to the library. Reasonable projects for first timers include:

• Filing out a bug report.
• Writing an example
• Improving the tests.
• Finding typos in the documentation.

By contributing code to scikit-fem, you are agreeing to release it under BSD-3-Clause, see LICENSE.md.

Citing the library

You may use the following BibTeX entry:

@article{skfem2020,
  doi = {10.21105/joss.02369},
  url = {https://doi.org/10.21105/joss.02369},
  year = {2020},
  publisher = {The Open Journal},
  volume = {5},
  number = {52},
  pages = {2369},
  author = {Tom Gustafsson and G. D. McBain},
  title = {scikit-fem: A Python package for finite element assembly},
  journal = {Journal of Open Source Software}
}

In literature

The library has been used in the preparation of the following scientific works:

• Gustafsson, T., Stenberg, R., & Videman, J. (2020). Nitsche's method for Kirchhoff plates. arXiv preprint arXiv:2007.00403.

• Gustafsson, T., & McBain, G. D. (2020). scikit-fem: A Python package for finite element assembly. Journal of Open Source Software, 52(5). Open access.

• Gustafsson, T., Stenberg, R., & Videman, J. (2020). On Nitsche's method for elastic contact problems. SIAM Journal on Scientific Computing, 42(2), B425–B446. arXiv preprint arXiv:1902.09312.

• Gustafsson, T., Stenberg, R., & Videman, J. (2019). Nitsche's Master-Slave Method for Elastic Contact Problems. arXiv:1912.08279.

• McBain, G. D., Mallinson, S. G., Brown, B. R., Gustafsson, T. (2019). Three ways to compute multiport inertance. The ANZIAM Journal, 60, C140–C155. Open access.

• Gustafsson, T., Stenberg, R., & Videman, J. (2019). Error analysis of Nitsche's mortar method. Numerische Mathematik, 142(4), 973–994. Open access.

• Gustafsson, T., Stenberg, R., & Videman, J. (2019). Nitsche's method for unilateral contact problems. Port. Math. 75, 189–204. arXiv preprint arXiv:1805.04283.

• Gustafsson, T., Stenberg, R. & Videman, J. (2018). A posteriori estimates for conforming Kirchhoff plate elements. SIAM Journal on Scientific Computing, 40(3), A1386–A1407. arXiv preprint arXiv:1707.08396.

• Gustafsson, T., Rajagopal, K. R., Stenberg, R., & Videman, J. (2018). An adaptive finite element method for the inequality-constrained Reynolds equation. Computer Methods in Applied Mechanics and Engineering, 336, 156–170. arXiv preprint arXiv:1711.04274.

• Gustafsson, T., Stenberg, R., & Videman, J. (2018). A stabilised finite element method for the plate obstacle problem. BIT Numerical Mathematics, 59(1), 97–124. arXiv preprint arXiv:1711.04166.

• Gustafsson, T., Stenberg, R., & Videman, J. (2017). Nitsche’s Method for the Obstacle Problem of Clamped Kirchhoff Plates. In European Conference on Numerical Mathematics and Advanced Applications, 407–415. Springer.

• Gustafsson, T., Stenberg, R., & Videman, J. (2017). A posteriori analysis of classical plate elements. Rakenteiden Mekaniikka, 50(3), 141–145. Open access.

Changelog

The format is based on Keep a Changelog, and this project adheres to Semantic Versioning.

Unreleased

Fixed

• Fix Mesh.validate for unsigned Mesh.t.

Added

• MeshTet.element_finder and MeshLine.element_finder for using InteriorBasis.interpolator.
• ElementTriCR, the nonconforming Crouzeix-Raviart element for Stokes flow.
• ElementTetCR, tetrahedral nonconforming Crouzeix-Raviart element.
• ElementTriHermite, an extension of ElementLineHermite to triangular meshes.

Deprecated

• L2_projection will be replaced by project.
• derivative will be replaced by project.

[2.1.1] - 2020-10-01

Fixed

• Further optimizations to Mesh3D.boundary_edges: tested to run on a laptop with over 10 million elements.

[2.1.0] - 2020-09-30

Fixed

• Mesh3D.boundary_edges (and, consequently, Basis.find_dofs) was slow and used lots of memory due to an exhaustive search of all edges

Added

• ElementHex2, a triquadratic hexahedral element
• MeshTri.init_circle, constructor for a circle mesh

[2.0.0] - 2020-08-21

Added

• Support for complex-valued forms: BilinearForm and LinearForm now take an optional argument dtype which defaults to np.float64 but can be also np.complex64
• Dofs.__or__ and Dofs.__add__, for merging degree-of-freedom sets (i.e. Dofs objects) using | and + operators
• Dofs.drop and Dofs.keep, for further filtering the degree-of-freedom sets

Removed

• Support for old-style decorators bilinear_form, linear_form, and functional (deprecated since 1.0.0)

Fixed

• FacetBasis did not initialize with ElementQuadP

Deprecated

• project will only support functions like lambda x: x[0] instead of lambda x, y, z: x in the future

[1.2.0] - 2020-07-07

Added

• Mesh.__add__, for merging meshes using + operator: duplicated nodes are joined
• ElementHexS2, a 20-node quadratic hexahedral serendipity element
• ElementLineMini, MINI-element for one-dimensional mesh

Fixed

• Mesh3D.boundary_edges was broken in case of hexahedral meshes
• skfem.utils.project did not work for ElementGlobal

Changed

• MeshQuad._splitquads aliased as MeshQuad.to_meshtri: should not be private

[1.1.0] - 2020-05-18

Added

• ElementTetMini, MINI-element for tetrahedral mesh

Fixed

• Mesh3D.boundary_edges incorrectly returned all edges where both nodes are on the boundary

[1.0.0] - 2020-04-22

Added

• New-style form constructors BilinearForm, LinearForm, and Functional
• skfem.io.json for serialization of meshes to/from json-files
• ElementLinePp, p-th order one-dimensional elements
• ElementQuadP, p-th order quadrilateral elements
• ElementQuadDG for transforming quadrilateral H1 elements to DG elements
• ElementQuadBFS, Bogner-Fox-Schmit element for biharmonic problems
• ElementTriMini, MINI-element for Stokes problems
• ElementComposite for using multiple elements in one bilinear form
• ElementQuadS2, quadratic Serendipity element
• ElementLineHermite, cubic Hermite element for Euler-Bernoulli beams
• Mesh.define_boundary for defining named boundaries
• Basis.find_dofs for finding degree-of-freedom indices
• Mesh.from_basis for defining high-order meshes
• Basis.split for splitting multicomponent solutions
• MortarMapping with basic support for mortar methods in 2D
• Basis constructors now accept quadrature keyword argument for specifying a custom quadrature rule

Deprecated

• Old-style form constructors bilinear_form, linear_form, and functional.

Changed

• Basis.interpolate returns DiscreteField objects instead of ndarray tuples
• Basis.interpolate works now properly for vectorial and high-order elements by interpolating all components and higher order derivatives
• Form.assemble accepts now any keyword arguments (with type DiscreteField) that are passed over to the forms
• Renamed skfem.importers to skfem.io
• Renamed skfem.models.helpers to skfem.helpers
• skfem.utils.solve will now expand also the solutions of eigenvalue problems