Single integral of the laplace kernel

Table of Contents

• Introduction
• The C++ code
  • Mesh generation
  • Function space and weighted integral
  • Results

Introduction

This tutorial explains how to use NiHu to compute a single singular integral of the form

\( \displaystyle I = \int_{S} K({\bf x}_0,{\bf y}) dS_y, \)

where \( S \) is a square surface

\( \displaystyle S = \left\{{\bf x} = (\xi,\eta,0) : 0 \le \xi \le 1, 0 \le \eta \le 1 \right\} \),

\( {\bf x}_0 \) is a point located at the surface's center \( (.5, .5, 0) \), and the integrand kernel \( K \) is the fundamental solution of the Laplace equation in 3D

\( \displaystyle K({\bf x},{\bf y}) = 1 / 4\pi r, \quad r = |{\bf x}-{\bf y}| . \)

The analytical result of the singular integral is

\( \displaystyle I = \log \left(1+\sqrt{2}\right) / \pi . \)

The C++ code

Mesh generation

The mesh of our surface domain will consist of only one 4-noded quadrilateral (NiHu::quad_1_elem) element:

    dMatrix nodes(4,3);
    nodes << 0., 0., 0., /*|*/ 1., 0., 0., /*|*/ 1., 1., 0., /*|*/ 0., 1., 0.;
    uMatrix elements(1, 1+4);
    elements << NiHu::quad_1_elem::id, 0, 1, 2, 3;
    auto msh = NiHu::create_mesh(nodes, elements, NiHu::quad_1_tag());

Function space and weighted integral

The single integral is generalised as a double weighted integral where one of the weighting functions is a Dirac delta function and the other is constant:

\( \displaystyle I = \int_{S} K({\bf x}_0,{\bf y}) dS_y = \int_{S} \delta({\bf x}_0) \int_{S} K({\bf x},{\bf y}) 1({\bf x}) dS_y dS_x = \left< \delta_{{\bf x}_0}, \left(\mathcal{K}1\right)_S \right>_S \)

The corresponding function spaces are instantiated as

    auto const &trial = NiHu::constant_view(msh);
    auto const &test = NiHu::dirac(trial);

As our mesh contains one element, and the piecewise constant function space contains one DOF per element, the resulting weighted integral matrix will contain one entry.

    dMatrix A(1, 1);
    A.setZero();
 
    auto K = NiHu::create_integral_operator(NiHu::laplace_3d_SLP_kernel());
    A << test * K[trial];

Results

The integration result is printed and compared to the analytical solution.

    std::cout << "WR matrix: " << A << std::endl;
 
    double anal = std::log(1. + root_two) / pi;
    std::cout << "log10 error = " << std::log10(std::abs(A.sum() / anal - 1.)) << std::endl;

The complete source of the tutorial is found here: tutorial/laplace_single_integral.cpp