Notebook

Exercise 2: Scattering from a pressure-release sphere¶

The problem¶

In the first and second tutorials, we looked at two formulations for a rigid scattering problem.

In this exercise, you will write your own code to solve a pressure-release scattering problem. As in the tutorials, we will use a unit sphere for $\Omega$, and we define the incident wave by

$$ p_{\text{inc}}(\mathbf x) = \mathrm{e}^{\mathrm{i} k x_0}. $$

where $\mathbf x = (x_0, x_1, x_2)$.

Acoustic waves are governed by the Helmholtz equation:

$$ \Delta p_\text{total} + k^2 p_\text{total} = 0, \quad \text{ in } \mathbb{R}^3 \backslash \Omega, $$

where $p_\text{total}$ is the total pressure. We can split $p_\text{total}$ into incident and scattered pressures by writing $p_\text{s}+p_\text{inc}$. The scattered pressure ($p_\text{s}$) satisfies the Sommerfeld radiation condition

$$ \frac{\partial p_\text{s}}{\partial r}-\mathrm{i}kp_\text{s}=o(r^{-1}) $$

when $r:=|\mathbf{x}|\rightarrow\infty$.

For our problem, we impose a Dirichlet boundary condition:

$$ p_\text{total}=0, \quad \text{ on } \Gamma, $$

where $\Gamma$ is the surface of the sphere $\Omega$.

The formulation¶

We use a formulation based on the direct formulation in the second tutorial.

Representation formula¶

For this problem, we use the following representation formula:

$$ p_\text{s} = \mathcal{D}u -\mathcal{S}\lambda, $$

where $\mathcal{S}$ is the single layer potential operator; $\mathcal{D}$ is the double layer potential operator; $u$ is the value (or trace) of $p_\text{s}$ on the surface $\Gamma$; and $\lambda$ is the normal derivative of $p_\text{s}$ on the surface $\Gamma$.

For this problem, our boundary condition tells us that $u=-p_\text{inc}$ on $\Gamma$.

Boundary integral equation¶

For this problem, we want to solve the following boundary integral equation:

$$ \mathsf{S}\lambda = -(\mathsf{D} - \tfrac{1}{2}\mathsf{I})p_\text{inc}, $$

where $\mathsf{S}$ is the single layer boundary operator; $\mathsf{D}$ is the double layer boundary operator, and $\mathsf{I}$ is the identity operator.

Solving with Bempp¶

Your task is to adapt and combine the two example codes in first and second tutorials to solve this problem and plot a slice of the solution at $z=0$

To get you started, I've copies the first few lines (which were the same in both examples) into the cell below.

In [1]:

%matplotlib inline
import bempp.api
from bempp.api.operators.boundary import helmholtz, sparse
from bempp.api.operators.potential import helmholtz as helmholtz_potential
from bempp.api.linalg import gmres
import numpy as np
from matplotlib import pyplot as plt

k = 8.

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

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

What next?¶

After attempting this exercises, you should read tutorial 3.