Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion .github/workflows/pythonapp.yml
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ jobs:
run: |
pip install mypy numpy
pip install git+https://github.com/FEniCS/ufl.git
mypy python/basix
mypy python/basix demo/python

build:
name: Build and test
Expand Down
9 changes: 8 additions & 1 deletion demo/python/demo_create_and_tabulate.py
Original file line number Diff line number Diff line change
Expand Up @@ -8,11 +8,18 @@
#
# First, we import Basix and Numpy.

import typing # For type checking

import numpy as np
import numpy.typing as npt

import basix
from basix import CellType, ElementFamily, LagrangeVariant

# Alias for type casting to maintain readability

FloatArray = npt.NDArray[np.float64]

# Next, we create a degree 4 Lagrange element on a quadrilateral using the function
# `create_element`. The first input is the element family: for Lagrange elements,
# we use `ElementFamily.P`. The second input is the cell type. The third
Expand All @@ -38,7 +45,7 @@
# functions (and no derivatives). We pass in the points as the second input.

points = np.array([[0.0, 0.0], [0.1, 0.1], [0.2, 0.3], [0.3, 0.6], [0.4, 1.0]])
tab = lagrange.tabulate(0, points)
tab = typing.cast(FloatArray, lagrange.tabulate(0, points))
print(tab)
print(tab.shape)

Expand Down
56 changes: 37 additions & 19 deletions demo/python/demo_custom_element.py
Original file line number Diff line number Diff line change
Expand Up @@ -7,11 +7,20 @@
#
# First, we import Basix and Numpy.

import typing # For type checking

import numpy as np
import numpy.typing as npt

import basix
from basix import CellType, LatticeType, MapType, PolynomialType, PolysetType, SobolevSpace

# Aliases for type casting to maintain readability

FloatArray = npt.NDArray[np.float64]
FloatingArray = npt.NDArray[np.floating]
QuadratureRule = tuple[FloatArray, FloatArray]

Comment on lines +20 to +23
Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

These are quite verbose changes to every demo for the sake of typing. What is the issue raised without them?

Copy link
Copy Markdown
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for the feedback! The demos had 305 errors with the mypy check which I've attached in a text file: demo_mypy_errors.txt

The majority of errors across all of the demos come from Basix functions returning broad ArrayLike types (these functions include make_quadrature, create_lattice, tabulate, tabulate_polynomials). Mypy is treating these values as large unions (buffers, strings, nested sequences, etc.), which leads to errors when the demos use them as NumPy arrays, for example:

  • slicing/indexing (pts[:, 0], poly[i, :])
  • arithmetic (f * poly * wts)
  • passing them to functions expecting NDArray

demo_custom_element.py had the most errors (199) which also included:

  1. Empty nested lists such as x = [[], [], [], []] and M = [[], [], [], []] are inferred as list[list[Never]], so mypy rejects later .append(np.array(...)) calls and the later create_custom_element(..., x, M, ...) calls
  2. There are also a few smaller API typing mismatches (using lists where tuples are expected).

The lines you commented on specifically are not necessary, I added them thinking it would keep things more readable but inline type casting would work just fine!

Copy link
Copy Markdown
Contributor

@schnellerhase schnellerhase May 1, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Right thanks. These are however exactly the errors we would like to resolve with this extended check, so ignoring them by casting to the expected types is not the desired resolution strategy. They show that the current type hints are insufficient. Rather we want to ensure that the type hints provided align with the use cases in the demos. For an example to resolve the errors related to tabulate_polynomials and make_quadrature see https://github.com/FEniCS/basix/tree/schnellerhase/fix-some-types (feel free to incorporate).

# Lagrange element with bubble
# ============================
#
Expand Down Expand Up @@ -76,11 +85,14 @@
# the largest degree that the integrand will be, so these integrals will
# be exact).

pts, wts = basix.make_quadrature(CellType.quadrilateral, 4)
poly = basix.tabulate_polynomials(PolynomialType.legendre, CellType.quadrilateral, 2, pts)
x = pts[:, 0]
y = pts[:, 1]
f = x * (1 - x) * y * (1 - y)
pts, wts = typing.cast(QuadratureRule, basix.make_quadrature(CellType.quadrilateral, 4))
poly = typing.cast(
FloatArray,
basix.tabulate_polynomials(PolynomialType.legendre, CellType.quadrilateral, 2, pts),
)
x_coord = pts[:, 0]
y_coord = pts[:, 1]
f = x_coord * (1 - x_coord) * y_coord * (1 - y_coord)
for i in range(9):
wcoeffs[4, i] = sum(f * poly[i, :] * wts)

Expand All @@ -100,7 +112,7 @@
#
# The shape of each of the point lists is (number of points, dimension).

x = [[], [], [], []]
x: list[list[FloatingArray]] = [[], [], [], []]
x[0].append(np.array([[0.0, 0.0]]))
x[0].append(np.array([[1.0, 0.0]]))
x[0].append(np.array([[0.0, 1.0]]))
Expand All @@ -121,7 +133,7 @@
# The shape of each matrix is (number of DOFs, value size, number of
# points, number of derivatives).

M = [[], [], [], []]
M: list[list[FloatingArray]] = [[], [], [], []]
for _ in range(4):
M[0].append(np.array([[[[1.0]]]]))
M[2].append(np.array([[[[1.0]]]]))
Expand Down Expand Up @@ -156,7 +168,7 @@

element = basix.create_custom_element(
CellType.quadrilateral,
[],
(),
Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is a correction we should make.

Related: it would be good to update the basix functions to accept any Sequence type and convert to tuples inside the functions. I've opened #1010 for this and labelled it good first issue

wcoeffs,
x,
M,
Expand Down Expand Up @@ -213,13 +225,16 @@
wcoeffs[0, 0] = 1
wcoeffs[1, 3] = 1

pts, wts = basix.make_quadrature(CellType.triangle, 2)
poly = basix.tabulate_polynomials(PolynomialType.legendre, CellType.triangle, 1, pts)
x = pts[:, 0]
y = pts[:, 1]
pts, wts = typing.cast(QuadratureRule, basix.make_quadrature(CellType.triangle, 2))
poly = typing.cast(
FloatArray,
basix.tabulate_polynomials(PolynomialType.legendre, CellType.triangle, 1, pts),
)
x_coord = pts[:, 0]
y_coord = pts[:, 1]
for i in range(3):
wcoeffs[2, i] = sum(x * poly[i, :] * wts)
wcoeffs[2, 3 + i] = sum(y * poly[i, :] * wts)
wcoeffs[2, i] = sum(x_coord * poly[i, :] * wts)
wcoeffs[2, 3 + i] = sum(y_coord * poly[i, :] * wts)

# Interpolation
# -------------
Expand All @@ -228,11 +243,11 @@
# the element are integrals. We begin by defining a degree 1 quadrature rule on an interval.
# This quadrature rule will be used to integrate on the edges of the triangle.

pts, wts = basix.make_quadrature(CellType.interval, 1)
pts, wts = typing.cast(QuadratureRule, basix.make_quadrature(CellType.interval, 1))

# The points associated with each edge are calculated by mapping the quadrature points to each edge.

x = [[], [], [], []]
x = typing.cast(list[list[FloatingArray]], [[], [], [], []])
for _ in range(3):
x[0].append(np.zeros((0, 2)))
x[1].append(np.array([[1 - p[0], p[0]] for p in pts]))
Expand All @@ -245,7 +260,7 @@
# edge, and no extra derivatives are used. The entries of these matrices are the quadrature weights
# multiplied by the normal directions.

M = [[], [], [], []]
M = typing.cast(list[list[FloatingArray]], [[], [], [], []])
for _ in range(3):
M[0].append(np.zeros((0, 2, 0, 1)))
for normal in [[-1, -1], [-1, 0], [0, 1]]:
Expand All @@ -260,7 +275,7 @@

element = basix.create_custom_element(
CellType.triangle,
[2],
(2,),
wcoeffs,
x,
M,
Expand All @@ -278,5 +293,8 @@

rt = basix.create_element(basix.ElementFamily.RT, CellType.triangle, 1)

points = basix.create_lattice(CellType.triangle, 1, LatticeType.equispaced, True)
points = typing.cast(
FloatArray,
basix.create_lattice(CellType.triangle, 1, LatticeType.equispaced, True),
)
assert np.allclose(rt.tabulate(0, points), element.tabulate(0, points))
21 changes: 17 additions & 4 deletions demo/python/demo_custom_element_conforming_cr.py
Original file line number Diff line number Diff line change
Expand Up @@ -7,13 +7,20 @@
#
# First, we import Basix and Numpy.

import typing # For type checking

import matplotlib.pyplot as plt
import numpy as np
import numpy.typing as npt
from mpl_toolkits import mplot3d # noqa: F401

import basix
from basix import CellType, LatticeType, MapType, PolynomialType, PolysetType, SobolevSpace

# Aliases for type casting to maintain readability

FloatArray = npt.NDArray[np.float64]

# Conforming CR element on a triangle
# ===================================
#
Expand Down Expand Up @@ -125,11 +132,14 @@ def create_ccr_triangle(degree):

# We now visualise the basis functions of the element we have created.

pts = basix.create_lattice(CellType.triangle, 30, LatticeType.equispaced, True)
pts = typing.cast(
FloatArray,
basix.create_lattice(CellType.triangle, 30, LatticeType.equispaced, True),
)

x = pts[:, 0]
y = pts[:, 1]
z = e.tabulate(0, pts)[0]
z = typing.cast(FloatArray, e.tabulate(0, pts))[0]

fig = plt.figure(figsize=(8, 8))

Expand All @@ -151,11 +161,14 @@ def create_ccr_triangle(degree):

e = create_ccr_triangle(3)

pts = basix.create_lattice(CellType.triangle, 30, LatticeType.equispaced, True)
pts = typing.cast(
FloatArray,
basix.create_lattice(CellType.triangle, 30, LatticeType.equispaced, True),
)

x = pts[:, 0]
y = pts[:, 1]
z = e.tabulate(0, pts)[0]
z = typing.cast(FloatArray, e.tabulate(0, pts))[0]

fig = plt.figure(figsize=(11, 8))

Expand Down
19 changes: 15 additions & 4 deletions demo/python/demo_dof_transformations.py
Original file line number Diff line number Diff line change
Expand Up @@ -22,11 +22,19 @@
#
# First, we import Basix and Numpy.

import typing # For type checking

import numpy as np
import numpy.typing as npt

import basix
from basix import CellType, ElementFamily, LagrangeVariant, LatticeType

# Aliases for type casting to maintain readability

FloatArray = npt.NDArray[np.float64]
IntArray = npt.NDArray[np.int_]

# Degree 5 Lagrange element
# =========================
#
Expand Down Expand Up @@ -129,7 +137,10 @@
# To demonstrate how these transformations can be used, we create a
# lattice of points where we will tabulate the element.

points = basix.create_lattice(CellType.tetrahedron, 5, LatticeType.equispaced, True)
points = typing.cast(
FloatArray,
basix.create_lattice(CellType.tetrahedron, 5, LatticeType.equispaced, True),
)

# If (for example) the direction of edge 2 in the physical cell does
# not match its direction on the reference, then we need to adjust the
Expand All @@ -145,10 +156,10 @@
# and over the value size. For each of these values, we apply the
# transformation matrix to the relevant DOFs.

data = nedelec.tabulate(0, points)
data = typing.cast(FloatArray, nedelec.tabulate(0, points))

transformation = nedelec.entity_transformations()["interval"][0]
dofs = nedelec.entity_dofs[1][2]
transformation = typing.cast(FloatArray, nedelec.entity_transformations()["interval"][0])
dofs = typing.cast(IntArray, np.asarray(nedelec.entity_dofs[1][2]))

for point in range(data.shape[1]):
for dim in range(data.shape[3]):
Expand Down
24 changes: 18 additions & 6 deletions demo/python/demo_facet_integral.py
Original file line number Diff line number Diff line change
Expand Up @@ -13,11 +13,20 @@
#
# We start by importing Basis and Numpy.

import typing # For type checking

import numpy as np
import numpy.typing as npt

import basix
from basix import CellType, ElementFamily, LagrangeVariant

# Aliases for type casting to maintain readability

FloatArray = npt.NDArray[np.float64]
IntArray = npt.NDArray[np.int_]
QuadratureRule = tuple[FloatArray, FloatArray]

# We define a degree 3 Lagrange space on a tetrahedron.

lagrange = basix.create_element(
Expand All @@ -28,7 +37,7 @@
# rule on a triangle. We use an order 3 rule so that we can integrate the
# basis functions in our space exactly.

points, weights = basix.make_quadrature(CellType.triangle, 3)
points, weights = typing.cast(QuadratureRule, basix.make_quadrature(CellType.triangle, 3))

# Next, we must map the quadrature points to our facet. We use the function
# `geometry` to get the coordinates of the vertices of the tetrahedron, and
Expand All @@ -39,8 +48,11 @@
#
# Using this information, we can map the quadrature points to the facet.

vertices = basix.geometry(CellType.tetrahedron)
facet = basix.cell.sub_entity_connectivity(CellType.tetrahedron)[2][0][0]
vertices = typing.cast(FloatArray, basix.geometry(CellType.tetrahedron))
facet = typing.cast(
IntArray,
basix.cell.sub_entity_connectivity(CellType.tetrahedron)[2][0][0],
)
mapped_points = np.array(
[
vertices[facet[0]] * (1 - x - y) + vertices[facet[1]] * x + vertices[facet[2]] * y
Expand All @@ -64,15 +76,15 @@
# We then multiply the three derivatives of the basis function by
# the three components of the normal.

normal = basix.cell.facet_outward_normals(CellType.tetrahedron)[0]
tab = lagrange.tabulate(1, mapped_points)[1:, :, 5, 0]
normal = typing.cast(FloatArray, basix.cell.facet_outward_normals(CellType.tetrahedron)[0])
tab = typing.cast(FloatArray, lagrange.tabulate(1, mapped_points))[1:, :, 5, 0]
normal_deriv = tab[0] * normal[0] + tab[1] * normal[1] + tab[2] * normal[2]

# As our facet is not the reference triangle, we must multiply the
# integrand by the norm of the Jacobian. We compute this by taking the
# cross product of the two columns of the Jacobian, and then compute
# the integral.

jacobian = basix.cell.facet_jacobians(CellType.tetrahedron)[0]
jacobian = typing.cast(FloatArray, basix.cell.facet_jacobians(CellType.tetrahedron)[0])
size_jacobian = np.linalg.norm(np.cross(jacobian[:, 0], jacobian[:, 1]))
print(np.sum(normal_deriv * weights) * size_jacobian)
20 changes: 16 additions & 4 deletions demo/python/demo_lagrange_variants.py
Original file line number Diff line number Diff line change
Expand Up @@ -25,11 +25,18 @@
#
# We begin by importing Basix and Numpy.

import typing # For type checking

import numpy as np
import numpy.typing as npt

import basix
from basix import CellType, ElementFamily, LagrangeVariant, LatticeType

# Alias for type casting

FloatArray = npt.NDArray[np.float64]

# In this demo, we consider Lagrange elements defined on a triangle. We start
# by creating a degree 15 Lagrange element that uses equally spaced points.
# This element will exhibit Runge's phenomenon, so we expect a large Lebesgue
Expand All @@ -55,8 +62,11 @@
#
# As expected, the value is large.

points = basix.create_lattice(CellType.triangle, 50, LatticeType.equispaced, True)
tab = lagrange.tabulate(0, points)[0]
points = typing.cast(
FloatArray,
basix.create_lattice(CellType.triangle, 50, LatticeType.equispaced, True),
)
tab = typing.cast(FloatArray, lagrange.tabulate(0, points))[0]
print(max(np.sum(np.abs(tab), axis=0)))

# A Lagrange element with a lower Lebesgue constant can be created by placing
Expand All @@ -70,12 +80,14 @@
# for the equally spaced element.

gll = basix.create_element(ElementFamily.P, CellType.triangle, 15, LagrangeVariant.gll_warped)
print(max(np.sum(np.abs(gll.tabulate(0, points)[0]), axis=0)))
gll_tab = typing.cast(FloatArray, gll.tabulate(0, points))[0]
print(max(np.sum(np.abs(gll_tab), axis=0)))

# An even lower Lebesgue constant can be obtained by placing the DOF points
# at GLL points mapped onto a triangle following the method proposed in
# `Recursive, Parameter-Free, Explicitly Defined Interpolation Nodes for
# Simplices (Isaac, 2020) <https://doi.org/10.1137/20M1321802>`_.

gll2 = basix.create_element(ElementFamily.P, CellType.triangle, 15, LagrangeVariant.gll_isaac)
print(max(np.sum(np.abs(gll2.tabulate(0, points)[0]), axis=0)))
gll2_tab = typing.cast(FloatArray, gll2.tabulate(0, points))[0]
print(max(np.sum(np.abs(gll2_tab), axis=0)))
Loading
Loading