This tutorial shows how to create and use a Calderón projectors and multitrace operators for Maxwell problems, as described in section 4 of Scroggs et al (2017).
For Maxwell problems, we define the multitrace operator $\mathsf{A}$ by
$$\mathsf{A}=\begin{bmatrix}\mathsf{H}&\mathsf{E}\\-\mathsf{E}&\mathsf{H}\end{bmatrix},$$where $\mathsf{H}$ and $\mathsf{E}$ are the magnetic and electric field boundary operators.
The interior Calderón projector $\mathcal{C}^\text{--}$ is defined by
$$\mathcal{C}^\text{--}=\tfrac12\mathsf{Id}+\mathsf{A},$$and the exterior Calderón projector $\mathcal{C}^\text{+}$ is defined by
$$\mathcal{C}^\text{+}=\tfrac12\mathsf{Id}-\mathsf{A}.$$If $\mathbf{E}$ and $\mathbf{H}$ are valid Cauchy data for a exterior Maxwell problem, then
$$\mathcal{C}^\text{+}\begin{bmatrix}\mathbf{E}\\\mathbf{H}\end{bmatrix}=\begin{bmatrix}\mathbf{E}\\\mathbf{H}\end{bmatrix}.$$It follows from the property $\left[\mathcal{C}^\text{+}\right]^2=\mathcal{C}^\text{+}$ that for any tangential functions $\mathbf{A}$ and $\mathbf{B}$, the functions $\mathbf{E}$ and $\mathbf{H}$ defined by
$$\begin{bmatrix}\mathbf{E}\\\mathbf{H}\end{bmatrix}=\mathcal{C}^\text{+}\begin{bmatrix}\mathbf{A}\\\mathbf{B}\end{bmatrix}$$are valid Cauchy data.
Let $\mathsf{B}$ and $\mathsf{C}$ be two boundary operators. If the product $\mathsf{BC}$ were discretised, the result would $B_hM^{-1}C_h$, where $M$ is the identity matrix between the range
and dual_to_range
spaces of $\mathsf{C}$.
If standard Raviart–Thomas (RT) and Nédélec (NC) basic functions are used to discretise $\mathsf{C}$, the matrix $M$ will be numerically singular. To avoid this problem, a stable pairing of RT and Buffa–Christiansen (BC) functions should be used.
In Bempp, the function bempp.api.operators.boundary.maxwell.multitrace_operator
will automatically use the correct spaces so that all terms in the product $\mathsf{A}^2$ are stable. This is discussed in more detail in section 4 of Scroggs et al (2017).
First, we import everything we need and set the wavenumber.
import bempp.api
from bempp.api.operators.boundary import maxwell
from bempp.api.operators.boundary import sparse
import numpy as np
k = 2
We make the meshed cube on which we will define our operators.
grid = bempp.api.shapes.cube(h=0.1)
Next, we make the multitrace operator and identity on the grid. Here, we must rell the multitrace identity to use vector "maxwell"
spaces instead of the default scalar spaces used for Laplace and Helmholtz.
multitrace = maxwell.multitrace_operator(grid, k)
identity = sparse.multitrace_identity(grid, spaces='maxwell')
We define the exterior Calderón projector, $\mathcal{C}^\text{+} = \tfrac12\mathsf{Id}-\mathsf{A}$.
calderon = 0.5 * identity - multitrace
We define grid functions using a Python function. These functions will not be valid Cauchy data.
def tangential_trace(x, n, domain_index, result):
result[:] = np.cross(np.array([1,0,0]),n)
electric_trace = bempp.api.GridFunction(
space=calderon.domain_spaces[0],
fun=tangential_trace,
dual_space=calderon.dual_to_range_spaces[0])
magnetic_trace = bempp.api.GridFunction(
space=calderon.domain_spaces[1],
fun=tangential_trace,
dual_space=calderon.dual_to_range_spaces[1])
We now apply the Calderón projector to the functions twice.
traces_1 = calderon * [electric_trace, magnetic_trace]
traces_2 = calderon * traces_1
As explained above, the functions in traces_1
should be valid Cauchy data. Additionally, traces_2
should be equal to traces_1
because $\left[\mathcal{C}^\text{+}\right]^2=\mathcal{C}^\text{+}$. To verify this, we plot the relative difference between the functions.
electric_error = (traces_2[0] - traces_1[0]).l2_norm() / traces_1[0].l2_norm()
magnetic_error = (traces_2[1] - traces_1[1]).l2_norm() / traces_1[1].l2_norm()
print("Electric error is ", electric_error)
print("Magnetic error is ", magnetic_error)
As expected, the errors are small.
traces_1[0].plot()