Welcome to modepy’s documentation!¶

modepy helps you create well-behaved high-order discretizations on simplices (i.e. triangles and tetrahedra) and tensor products of simplices (i.e. squares, cubes, prisms, etc.). These are a key building block for high-order unstructured discretizations, as often used in a finite element context.

The basic objects that modepy manipulates are functions on a shape (or reference domain). For example, it supplies an orthonormal basis on triangles (shown below).

Proriol-Koornwinder-Dubiner orthogonal (PKDO) basis functions of order 3

The file that created this plot is included in the modepy distribution as examples/plot-basis.py.

Its roots closely followed the approach taken in the Hesthaven and Warburton book [HesthavenWarburton2007], but much has been added beyond that basic functionality.

[HesthavenWarburton2007]

J. S. Hesthaven and T. Warburton (2007). Nodal Discontinuous Galerkin Methods: Algorithms, Analysis, and Applications (1st ed.). doi:10.1007/978-0-387-72067-8. (source code)

Example¶

Here’s an idea of code that uses modepy:

Example code that constructs a prism domain and computes a weak derivative on it.¶
 3import numpy as np
 4
 5import modepy as mp
 6
 7
 8# Define the shape on which we will operate
 9line = mp.Simplex(1)
10triangle = mp.Simplex(2)
11prism = mp.TensorProductShape((triangle, line))
12
13assert prism.dim == 3
14
15# Define a function space for the prism
16n = 12
17space = mp.TensorProductSpace((mp.PN(triangle.dim, n), mp.PN(line.dim, n)))
18
19assert space.order == n
20assert space.spatial_dim == 3
21assert space.space_dim == (n + 1) * (n + 1) * (n + 2) // 2
22
23# Define a basis function for the prism
24basis = mp.orthonormal_basis_for_space(space, prism)
25
26# Define a point set for the prism
27nodes = mp.edge_clustered_nodes_for_space(space, prism)
28
29# Define a quadrature rule for the prism
30quadrature = mp.TensorProductQuadrature([
31    mp.VioreanuRokhlinSimplexQuadrature(n, triangle.dim),
32    mp.LegendreGaussQuadrature(n),
33])
34
35# Define a bilinear form: weak derivative in the x direction
36i = 0
37weak_d = mp.nodal_quadrature_bilinear_form_matrix(
38    quadrature,
39    test_functions=basis.derivatives(i),
40    trial_functions=basis.functions,
41    nodal_interp_functions_test=basis.functions,
42    nodal_interp_functions_trial=basis.functions,
43    input_nodes=nodes,
44    output_nodes=nodes,
45)
46inv_mass = mp.inverse_mass_matrix(basis, nodes)
47
48# Compute derivative
49f = 1.0 - np.sin(nodes[i]) ** 3
50f_ref = -3.0 * np.cos(nodes[i]) * np.sin(nodes[i]) ** 2
51f_approx = inv_mass @ weak_d.T @ f
52
53error = np.linalg.norm(f_approx - f_ref) / np.linalg.norm(f_ref)
54assert error < 1.0e-4

This file is included in the modepy distribution as examples/prism-forms.py.

modepy around the web¶

Contents¶