# Boundary elements method — calculation of solid angle

I am developing a BEM code based on a deal.ii tutorial. Consider the Poisson equation

$$\Delta u=-f\,,$$ and its Green's function $$G\left(\mathbf{x},\mathbf{x}'\right)$$ with the property $$\Delta G=-\delta\left(\mathbf{x}-\mathbf{x}'\right)\,.$$ Applying Green's 2nd identity leads to the well known boundary integral equation (BIE) $$\int\limits_{\Omega}u\left(\mathbf{x}\right)\delta\left(\mathbf{x}-\mathbf{x}'\right)\mathrm{d}V =\int\limits_{\Omega}f\left(\mathbf{x}\right)G\left(\mathbf{x},\mathbf{x}'\right)\mathrm{d}V +\oint\limits_{\partial\Omega}\mathbf{n}\cdot\left[\cdots\right]\mathrm{d}\Gamma$$ I am considering the integral on left hand side of the equation. It is well known that $$\int\limits_{\Omega}u\left(\mathbf{x}\right)\delta\left(\mathbf{x}-\mathbf{x}'\right)\mathrm{d}V=\alpha\left(\mathbf{x}'\right)u\left(\mathbf{x}'\right)\,,$$ where as the factor $$\alpha$$, i.e. the fraction of the solid angle, depends on the location of the point $$\mathbf{x}'$$. Typically it is given by $$\alpha\left(\mathbf{x}\right)=\begin{cases} 0 & \mathbf{x}\notin\Omega\\ 1 & \mathbf{x}\in\Omega\\ \frac{1}{2} & \mathbf{x}\in\Gamma & \text{smooth boundary}\\ \frac{\theta}{4\pi} & \mathbf{x}\in\Gamma & \text{3D corner} \end{cases}$$ You get these expressions by simple analytic evaluations of the integral. If the BIE is discretized, these analytic expressions may be not be appropiate any more. Therefore, the author of the tutorial states that the fraction of solid angle may be calculate by $$\alpha\left(\mathbf{x}'\right)=1+\oint\mathbf{n}\cdot\nabla G\left(\mathbf{x},\mathbf{x}'\right)\mathrm{d}\Gamma\,,$$ and in fact uses this term in the discretization.

Finally my question: How can be the equation above can be derived?

• I solved this. Tomorrow I will type the solution. – sebastian_g Jan 14 '16 at 18:15
• This question would best have been answered on the deal.II mailing list. In fact, Luca (the author of the program) may have answered it :-) – Wolfgang Bangerth Jan 14 '16 at 19:21
• But since you already have the answer, may I suggest you simply write it up as a comment in the file examples/step-34/doc/intro.dox so that others who read it have the benefit of having the information there? – Wolfgang Bangerth Jan 14 '16 at 19:22
• @WolfgangBangerth I cloned deal.ii's git repository and modified the file examples/step-34/doc/intro.dox. I am allowed to make a commit and push it to the online repository? – sebastian_g Jan 19 '16 at 10:24
• No, but you can extract the patch and send it to me/us. See for example git help format-patch to see how to do this. (In the future, the way to do it is to "fork" the repository on your own github account, "clone" your fork onto your own harddrive, make modifications there, push to your fork, and create a pull request from your fork to the main repository.) – Wolfgang Bangerth Jan 19 '16 at 14:02

If a homogeneous Neumann problem is considered, i.e. $f=0$ in $\Omega$ and $\mathbf{n}\cdot\nabla\phi$ on $\Gamma$, one solution is given by the constant function
$$\phi\left(\mathbf{x}\right)=\phi_0\,.$$
$$\int\limits_{\Omega}\phi\left(\mathbf{x}' \right)\delta\left(\mathbf{x}'-\mathbf{x}\right)\mathrm{d}V' =\phi_0\int\limits_{\Omega}\delta\left(\mathbf{x}'-\mathbf{x}\right)\mathrm{d}V' =-\phi_0\oint\limits_{\Gamma}\tilde{\mathbf{n}}\cdot \left[\nabla' G\left(\mathbf{x}'-\mathbf{x}\right)\right]\mathrm{d}\Gamma'\,,$$
$\tilde{\mathbf{n}}$ is the outward normal of the considered domain. Let the considered domain be the whole space substracted with an interior domain $\Gamma_0$. Then the boundary integral reads $$-\oint\limits_{\Gamma}\tilde{\mathbf{n}}\cdot \left[\nabla' G\left(\mathbf{x}'-\mathbf{x}\right)\right]\mathrm{d}\Gamma' =-\oint\limits_{\Gamma_\infty}\tilde{\mathbf{n}}\cdot \left[\nabla' G\left(\mathbf{x}'-\mathbf{x}\right)\right]\mathrm{d}\Gamma'+\\ +\oint\limits_{\Gamma_0}\mathbf{n}\cdot \left[\nabla' G\left(\mathbf{x}'-\mathbf{x}\right)\right]\mathrm{d}\Gamma'\,,$$ where $\mathbf{n}$ is the inward normal. It can be shown that the integral on $\Gamma_\infty$ is unity (use spherical coordinates). Using this result and the 2nd equation, we have $$\alpha\left(\mathbf{x}\right) =1+\oint\limits_{\Gamma_0}\mathbf{n}\cdot \left[\nabla' G\left(\mathbf{x}'-\mathbf{x}\right)\right]\mathrm{d}\Gamma'\,.$$