Function spaces

Similarly to finite element codes, function spaces form a central part of Bempp. To initialize a function space all we need is a grid object.

import bempp.api
grid = bempp.api.shapes.regular_sphere(3)

Let us now create a space of piecewise constant functions on the elements. This is the standard low-order space to represent Neumann data (that is normal derivatives) in boundary element methods.

space = bempp.api.function_space(grid, "DP", 0)

The first parameter of the function_space function is always the grid. The second gives the type of space, in this case "DP" for discontinuous polynomial, and the third parameter is the order of the space (0 for piecewise constant). We can now query the degrees of freedom (dofs) of the space.

number_of_global_dofs = space.global_dof_count
print("Number of global dofs: {0}".format(number_of_global_dofs))
Number of global dofs: 512

For this space, the number of dofs is equal to the number of elements on the mesh. This is not necessarily always the case. Let us now create a space of continuous, piecewise linear functions. This is the standard low-order space to represent Dirichlet data.

space = bempp.api.function_space(grid, "P", 1)
number_of_global_dofs = space.global_dof_count
print("Number of global dofs: {0}".format(number_of_global_dofs))
Number of global dofs: 258

The number of global dofs is now equal to the number of vertices in the grid.

Space types

Bempp supports the following spaces.

Scalar spaces

The following scalar spaces are well suited to Laplace and Helmholtz problems.

Identifier Space Description
"DP" Discontinuous polynomial spaces These are spaces of functions that are polynomial on each element but discontinuous across elements. The maximum order is 10.
"P" Polynomial spaces These are spaces of functions that are polynomial on each element and continuous across elements. The minimum order is 1. The maximum order is 10.
"DUAL" Piecewise polynomials on the dual mesh These spaces are implemented for order 0 only. These spaces form a stable dual pairing together with continuous, piecewise linear functions and are needed for certain opposite order preconditioners. When using these spaces, "B-DP" and "B-P" spaces (versions of "DP" and "P" defined on the barycentrically refined grid) should be used.

Vector spaces

The following vector spaces are well suited to Maxwell problems.

Identifier Space Description
"RT" Raviart–Thomas spaces These are spaces of div-conforming Raviart-Thomas basis functions. Currently, only order 0 Raviart-Thomas spaces are supported.
"NC" Nédélec spaces These are spaces of curl-conforming Nédélec basis functions. Currently, only order 0 Nédélec spaces are supported.
"BC" Buffa–Christansen dual spaces These are spaces of div-conforming and quasi-curl-conforming Buffa–Christiansen basis functions. There are heavily used to implement Calderón preconditioning. Currently, only order 0 Buffa–Christansen spaces are supported. As with "DUAL"spaces, the barycentric versions of spaces ("B-RT" and "B-NC") should be used in conjunction with these spaces.
"RBC" Rotated Buffa–Christansen dual spaces These are spaces of curl-conforming and quasi-div-conforming basis functions given by taking the cross product of "BC" basis functions with the normal to the surface.

Local and global dofs

An important concept for spaces are global and local degrees of freedom (dofs). Global dofs are the actual dofs of the space associated with the basis functions, while local dofs are the coefficients of the restrictions of the basis functions to individual elements. Consider for example a space of continuous, piecewise linear basis functions. Each global dof is associated with a vertex. The corresponding basis function lives on all elements that are adjacent to this vertex.

Sometimes it is useful to find out which global dofs the basis functions on an element contribute. This is shown in the following example.

space = bempp.api.function_space(grid, "P", 1)
elements = list(grid.leaf_view.entity_iterator(0))
element = elements[0]
global_dofs, weights = space.get_global_dofs(element, dof_weights=True)
print("Map from local to global dofs on element: {0}".format(global_dofs))
print("Dof weights: {0} ".format(weights))
Map from local to global dofs on element: [2, 66, 68]
Dof weights: [1.0, 1.0, 1.0]

The map from local to global dofs on the element denotes which local basis function is associated with which global basis function. In this example, the element has three local basis functions (associated with the three vertices of the element), and the values of global_dofs are the indices of the global basis functions that they map to. The array weights returns scalar multipliers that are the weights with which the local basis function contributes to the global basis function. Hence, to obtain the local coefficient of a local basis function the corresponding global dof coefficient needs to be multiplied with the weight. For most basis types the weight is always 1. But for example, for Raviart–Thomas basis functions it can differ. Weights are only returned if the parameter dof_weights=True is set.