Available operators

In this tutorial, we describe all operators available in Bempp. For details of the functionality of these operators see the operator tutorial.

Scalar non-local operators

Let \displaystyle g(\mathbf{x},\mathbf{y}) be a given Green’s function. We define the single layer (\mathcal{V}) and double layer (\mathcal{K}) potential operators as follows.

[\mathcal{V}\phi](\mathbf{x}) = \int_{\Gamma}g(\mathbf{x},\mathbf{y})\phi(\mathbf{y})\mathrm{d}s(\mathbf{y}),~\mathbf{x}\in\mathbb{R}^3\backslash\Gamma\nonumber\\[2mm]   [\mathcal{K}\phi](\mathbf{x}) = \int_{\Gamma}\frac{\partial g(\mathbf{x},\mathbf{y})}{\partial\nu(\mathbf{y})}\phi(\mathbf{y})\mathrm{d}s(\mathbf{y}),~\mathbf{x}\in\mathbb{R}^3\backslash\Gamma

From this we derive the following boundary operators:

Operator Formula
Single layer boundary operator \displaystyle [\mathsf{V}\phi](\mathbf{x}) = \int_{\Gamma}g(\mathbf{x},\mathbf{y})\phi(\mathbf{y})\mathrm{d}\mathbf{y},~\mathbf{x}\in\Gamma
Double layer boundary operator \displaystyle [\mathsf{K}\phi](\mathbf{x}) = \int_{\Gamma}\frac{\partial g(\mathbf{x},\mathbf{y})}{\partial\nu(\mathbf{y})}\phi(\mathbf{y})\mathrm{d}\mathbf{y},~\mathbf{x}\in\Gamma
Adjoint double layer boundary operator \displaystyle [\mathsf{K}'\phi](\mathbf{x}) = \int_{\Gamma}\frac{\partial g(\mathbf{x},\mathbf{y})}{\partial\nu(\mathbf{x})}\phi(\mathbf{y})\mathrm{d}\mathbf{y},~\mathbf{x}\in\Gamma
Hypersingular boundary operator \displaystyle [\mathsf{D}\phi](\mathbf{x}) = -\frac{\partial}{\partial \nu(\mathbf{x})}\int_{\Gamma}\frac{\partial g(\mathbf{x},\mathbf{y})}{\partial\nu(\mathbf{y})}\phi(\mathbf{y})\mathrm{d}\mathbf{y},~\mathbf{x}\in\Gamma

The actual implementation of boundary operators in Bempp is based on variational formulations. Foe example, the implementation of the single-layer boundary operator is given by

s(\phi, \psi) = \int_{\Gamma}\overline{\psi(\mathbf{x})}\int_{\Gamma}g(\mathbf{x},\mathbf{y})\phi(\mathbf{y})\mathrm{d}\mathbf{y}\mathrm{d}\mathbf{x}.

For the hypersingular operator, a slightly different formulation based on integration by parts of this variational form is used. Details can be found in the book Numerical Approximation Methods for Elliptic Boundary Value Problems by Steinbach.

The different types of boundary operators are dependent on the choice of the Green’s function \displaystyle g. The Green’s functions used in Bempp, plus the modules in which the boundary operators are defined, are given below.

PDE Green’s function Module
Laplace
\displaystyle -\Delta u = 0
g(\mathbf{x},\mathbf{y}) = \frac{1}{4\pi |\mathbf{x}-\mathbf{y}|} bempp.api.operators.boundary.laplace
Modified Helmholtz
\displaystyle -\Delta u+\omega^2 u = 0
g(\mathbf{x},\mathbf{y}) = \frac{\mathrm{e}^{-\omega |\mathbf{x}-\mathbf{y}|}}{4\pi|\mathbf{x}-\mathbf{y}|} bempp.api.operators.boundary.modified_helmholtz
Helmholtz
\displaystyle -\Delta u - k^2 u = 0
g(\mathbf{x},\mathbf{y}) = \frac{\mathrm{e}^{\mathrm{i}k |\mathbf{x}-\mathbf{y}|}}{4\pi|\mathbf{x}-\mathbf{y}|} bempp.api.operators.boundary.helmholtz

In each of the modules, the names of the corresponding functions are single_layer, double_layer, adjoint_double_layer and hypersingular.

The packages for the corresponding potential operators are:

PDE Module
Laplace bempp.api.operators.potential.laplace
Helmholtz bempp.api.operators.potential.helmholtz
Modified Helmholtz bempp.api.operators.potential.modified_helmholtz

The names of the corresponding potentials in each module are single_layer and double_layer.

For example, a Helmholtz hypersingular boundary operator is instantiated as

hyp = bempp.api.operators.boundary.helmholtz.hypersingular(
domain, range, dual_to_range, k)

Here, k is the wavenumber that appears in the PDE and Green’s function in the table above.

In several applications, in particular transmission problems it is useful to define the multitrace operator

\mathsf{A} = \begin{bmatrix} -\mathsf{K} & \mathsf{V}\\ \mathsf{D} & \mathsf{K}'\end{bmatrix}.

This operator has the property that

\left[\tfrac12\mathsf{Id} + \mathsf{A}\right]\begin{bmatrix}u \\ u_n \end{bmatrix} = \begin{bmatrix}u \\ u_n\end{bmatrix}

if \displaystyle u and \displaystyle u_n are the Dirichlet and Neumann data (respectively) of an interior solution to the corresponding PDE. The operator \displaystyle \tfrac12\mathsf{Id} + \mathsf{A} is called the interior Calderón projector. From this property, it follows directly that \displaystyle (2\mathsf{A})^2 = \mathsf{Id} and therefore \displaystyle \mathsf{A} is self-regularizing.

The exterior Calderón projector is defined as \displaystyle \tfrac12\mathsf{Id} - \mathsf{A}.

The following code snippet defines a multitrace operator for a Helmholtz problem with wavenumber \displaystyle k=1 and forms the interior Calderón projector.

grid = bempp.api.shapes.regular_sphere(3)
A = bempp.api.operators.helmholtz.multitrace_operator(grid, 1)
ident = bempp.api.operators.boundary.sparse.multitrace_identity(grid)
calderon = .5 * ident + A

The assembly of the multitrace operator \displaystyle \mathsf{A} requires the assembly of the single-layer and double-layer boundary operators using discontinuous piecewise linear basis functions on a barycentric refinement of the grid. This barycentric refinement has six times the number of elements as the original grid. This is necessary to assemble a basis for piecewise constant basis functions on the dual grid.

Scalar, sparse operators

Bempp implements two sparse operators acting on scalar spaces: the identity operator \displaystyle \mathsf{Id} and the Laplace-Beltrami operator \displaystyle -\Delta_{\Gamma}. They are given in their variational form as follows.

Operator Variational form Function
Identity Operator m(\phi, \psi) = \int_{\Gamma} \overline{\psi(\mathbf{y})} \phi(\mathbf{y})\mathrm{d}\mathbf{y} bempp.api.boundary.sparse.identity
Laplace-Beltrami Operator k(\phi, \psi) = \int_{\Gamma} \overline{\nabla_{\Gamma}\psi(\mathbf{y})}\cdot \nabla_{\Gamma}\phi(\mathbf{y})\mathrm{d}\mathbf{y} bempp.api.boundary.sparse.laplace_beltrami

Maxwell operators

Bempp supports the solution of the time-harmonic Maxwell equation of the form

\nabla\times\nabla\times \mathbf{E} - k^2 \mathbf{E} = 0.

We define the electric field (\mathcal{E}) and magnetic field (\mathcal{H}) potential operators by:

[\mathcal{E}\Phi](\mathbf{x}) = \mathrm{i}k\int_{\Gamma}g(\mathbf{x},\mathbf{y})\Phi(\mathbf{y})\mathrm{d}\mathbf{y}-\frac{1}{\mathrm{i}k}\nabla_\mathbf{x}\int_{\Gamma}g(\mathbf{x},\mathbf{y})(\nabla_{\Gamma}\cdot\Phi)(\mathbf{y})\mathrm{d}\mathbf{y}
[\mathcal{H}\Phi](\mathbf{x}) = \nabla_\mathbf{x}\times\int_{\Gamma}g(\mathbf{x},\mathbf{y})\phi(\mathbf{y})\mathrm{d}\mathbf{y}

The corresponding functions are bempp.api.operators.potential.maxwell.electric_field and bempp.api.operators.potential.maxwell.magnetic_field.

The definition of the electric field operator given above includes a factor \displaystyle \mathrm{i} in the numerator. This is less common, but has implementational advantages in Bempp.

Based on the electric and magnetic field potential operators we can derive the corresponding boundary operators. These are given in the following table.

Operator Varational form & function
Electric field boundary operator s(\Phi, \Psi\times\nu) = \int_{\Gamma}\int_{\Gamma}g(\mathbf{x},\mathbf{y})\left[-\mathrm{i}k\overline{\Psi(\mathbf{x})}\cdot\Phi(\mathbf{y})-\frac{1}{\mathrm{i}k}\left(\nabla_{\Gamma}\cdot\overline{\Psi}\right)(\mathbf{x})\left(\nabla_{\Gamma}\cdot\Phi\right)(\mathbf{y})\right]\mathrm{d}\mathbf{x}\mathrm{d}\mathbf{y}
bempp.operators.boundary.maxwell.electric_field
Magnetic field boundary operator c(\Phi, \Psi\times\nu) = \int_{\Gamma}\int_{\Gamma}\nabla_\mathbf{x}g(\mathbf{x},\mathbf{y})\left[\overline{\Psi(\mathbf{x})}\times \Phi(\mathbf{y})\right]\mathrm{d}\mathbf{x}\mathrm{d}\mathbf{y}
bempp.operators.boundary.maxwell.magnetic_field

In these variational forms, the domain functions \Phi and the range functions should belong div-conforming spaces, while the dual_to_range functions \Psi\times\nu should be curl-conforming. More details of these spaces can be found in the spaces tutorial and the Maxwell tutorials.

The Maxwell multitrace operator is defined by

\mathsf{A}=\begin{bmatrix}\mathsf{H}&\mathsf{E}\\-\mathsf{E}&\mathsf{H}\end{bmatrix}

and the interior and exterior Calderón projectors are defined as \displaystyle \tfrac12\mathsf{Id} + \mathsf{A} and \displaystyle \tfrac12\mathsf{Id} - \mathsf{A}. The Maxwell multitrace operator (with wavenumber 1) and interior Calderón projector may be created in Bempp using the following code.

grid = bempp.api.shapes.regular_sphere(3)
A = bempp.api.operators.maxwell.multitrace_operator(grid, 1)
ident = bempp.api.operators.boundary.sparse.multitrace_identity(grid,
use_maxwell_space=True)
calderon = .5 * ident + A

These operators in the multitrace operator are automatically discretised using stable pairings of Raviart–Thomas and Buffa–Christiansen basis functions. More details on this can be found in this tutorial.

Far Field Operators

Bempp implements far field operators for Maxwell and Helmholtz problems. These can be found in the modules bempp.api.operators.far_field.maxwell and bempp.api.operators.far_field.helmholtz.

Approximate DtN and NtD operators

Bempp implements the OSRC approximations to the Dirichlet-to-Neumann (DtN) and Neumann-to-Dirichlet (NtD) maps (see Antoine and Darbas (2007)). The corresponding functions are bempp.api.boundary.helmholtz.osrc_ntd and bempp.api.boundary.helmholtz.osrc_dtn.