%               USING Texmaker

documentclass{article}
usepackage{multicol}
usepackage[margin=0.5in]{geometry}
usepackage[T1]{fontenc}
usepackage{sectsty}
usepackage{amsmath}
usepackage{stfloats}
usepackage{graphicx}
usepackage{wrapfig}

usepackage{array}
usepackage{calc}
usepackage{mwe}

usepackage{polski}
usepackage[utf8]{inputenc}


sectionfont{fontsize{10}{10}selectfont}
subsectionfont{fontsize{10}{10}selectfont}
renewcommand{figurename}{Rys.}

begin{document}

begin{center}
  setlength{tabcolsep}{0pt}
  begin{tabular}{>{raggedleft}m{1cm}>{centering}m{dimexprtextwidth - 5cmrelax}>{raggedright}m{5cm}}
  includegraphics[width=2cm]{./top_left.png}%
  &&
  includegraphics[width=5cm]{./top_right_w.png}
  end{tabular}
  end{center}

begin{center}
Engineering Analysis with Boundary Elements 30 (2006) $371-381$ par
%setlength{parindent}{10ex}
begin{Large}
setlength{parindent}{24cm}
textbf{   The fast multipole boundary element method for \ potential problems: A tutorial}
end{Large}
end{center}

begin{center}
begin{Large}
Y.J. Liu $mathrm{a},*$, N. Nishimura $mathrm{b},1$
end{Large}

begin{small}
{it a Department of Mechanical, Industrial and Nuclear Engineering, University of Cincinnati, P. O. Box 210072, Cincinnati}, $OH45221-0072$, {it USA} $mathrm{b}$ {it Academic Center for Computing and Media Studies, Kyoto University, Kyoto 606-8501, Japan}

Received 6 May 2005; accepted 23 November 2005 Available online 15 March 2006
end{small}
end{center}


begin{flushleft}
 ------------------------------------------------------------------------------------------------------------------------------------------------------------------
end{flushleft}

section{Abstract}

indent indent The fast multipole method (FMM) has been regarded as one of the top 10 algorithms in scientific computing that were developed in the 20th century. Combined with the FMM, the boundary element method (BEM) can now solve large-scale problems with several million degrees of freedom on a desktop computer within hours. This opened up a wide range of applications for the BEM that has been hindered for many years by the lack of eficiencies in the solution process, although it has been regarded as superb in the modeling stage. However, understanding the fast multipole BEM is even more difficult as compared with the conventional BEM, because of the added complexities and different approaches in both FMM formulations and implementations. This paper is an introduction to the fast multipole BEM for potential problems, which is aimed to overcome this hurdle for people who are familiar with the conventional BEM and want to learn and adopt the fast multipole approach. The basic concept and main procedures in the FMM for solving boundary integral equations are described in detail using the $2mathrm{D}$ potential problem as an example. The structure of a fast multipole BEM program is presented and the source code is also made available that can help the development of fast multipole BEM codes for solving other problems. Numerical examples are presented to further demonstrate the eficiency, accuracy and potentials of the fast multipole BEM for solving large-scale problems.
begin{flushleft}
copyright 2006 Elsevier Ltd. All rights reserved.

{it Keywords}: Fast multipole method; Boundary element method; $2mathrm{D}$ potential problems
end{flushleft}

begin{flushleft}
 ------------------------------------------------------------------------------------------------------------------------------------------------------------------
end{flushleft}

begin{multicols}{2}
raggedcolumns
section{Introduction}
The boundary integral equation (BIE) formulations and their numerical solutions using boundary element method (BEM) for mechanics problems originated about 40 years ago. The 2D potential problem was first formulated in terms of a direct BIE and solved by Jaswon [1]. This work was later extended by Rizzo to the vector case—2D elastostatic problem [2]. Following these early works, extensive research efforts have been made for the development of the BIE/BEM. Some of the important textbooks and research volumes in English can be found in Refs. [3–9].
Although the BEM has enjoyed the reputation of easy meshing in modeling for many problems, its efficiency in solutions has been a serious problem for analyzing large-scale models that have emerged in applications with the availability of more powerful computers. For example, while the finite element method (FEM) has been used routinely to solve models with several millions of degrees of freedom (DOF’s), the BEM has been limited to solving problems with a few thousands DOF’s for many years.This is because the conventional BEM in general produces dense and non-symmetric matrices that, although smaller in sizes, requires O(N2) operations to compute the coefficients and another O(N3) operations to solve the system using direct solvers (where Nis the number of equations of the linear system).
In the mid of 1980s, Rokhlin and Greengard [10–12] pioneered the innovative fast multipole method (FMM) that can be used to accelerate the solutions of boundary integral equations by several folds, promising to reduce the CPU time in FMM accelerated BEM to O(N). However, it took almost a decade for the mechanics community to realize the potential of the FMM for the BEM. Some of the early research on FMM– BEM in applied mechanics can be found in Refs. [13–17], which show great promises of the FMM–BEM for solving large-scale engineering problems. Most recently, composite material models containing tens of thousands of fibers [18,19] and models of electromagnetic wave scatterings from a full aircraft at giga hertz frequencies [20] all have been solved successfully by using the FMM–BEM within hours and with moderate computing resources. A comprehensive review of the FMM BIE/BEM can be found in Ref. [21].
However, the use of the fast multipole method has increased the complexity in implementations of the BIE/BEM. Although there are many research papers dealing with the various subjects of the FMM–BEM, very few can be used as introductory materials to learn the new approach for researchers who are familiar with the conventional BEM or for new comers who are often frustrated with even the conventional BIE/BEM formulations and implementations. The research papers are often concise, without detailed information, which made the understanding of the FMM– BEM difficult. To fill the gap and hence promote the FMM– BEM in the BIE/BEM research community, we will provide in this paper an introduction to the FMM–BEM with complete information on formulations, discretizations and implemen-tations using the 2D potential problem as an example.
This paper is organized as follows: in Section 2, we review the BIE formulation for potential problems and the conven-tional BEM in discretization of the BIE. In Section 3, we present the complete formulations and algorithms used in the FMM–BEM for 2D potential problem. In Section 4, we discuss the programming issues in FMM–BEM using a Fortran code that is available from the authors with this paper. In Section 5, we show some example problems solved by using the provided code to demonstrate the efficiencies of the FMM–BEM for large-scale problems. We conclude the article with some discussions in Section 6.

section{BIE formulation and the conventional BEM approach}
The boundary integral equation formulation and its discretization using the conventional BEM for potential problems is summarized in this section. Consider the following Laplace equation governing potential problems (for example, a steady-stateheatconductionproblem)ina2DdomainV(Fig.1) \
begin {equation}
nabla^{2} phi(mathbf{x})=0, quad forall mathbf{x} in V
 end {equation}
begin{wrapfigure}{l}{linewidth}
vspace{-20pt}
begin{center}
includegraphics[width=linewidth]{./fig_1.png}
end{center}
vspace{-20pt}
caption{Figure}
end{wrapfigure}

indent \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \
under the boundary conditions \
begin {equation}
phi(mathbf{x})=overline{phi}(mathbf{x}), quad forall mathbf{x} in S_{1}
end {equation}
begin {equation}
q(mathbf{x}) equiv frac{partial phi}{partial n}(mathbf{x})=overline{q}(mathbf{x}), quad forall mathbf{x} in S_{2}
end {equation}
where f is the potential field in domain V, SZS1gS2 the boundary of V, n the outward normal, and the barred quantities indicate given values on the boundary.
The solution to the boundary value problem described by Eqs. (1)–(3) can be written as the following representation integral (see, e.g. Refs. [7–9])\
begin {equation}
phi(mathbf{x})=int_{S}[G(mathbf{x}, mathbf{y}) q(mathbf{y})-F(mathbf{x}, mathbf{y}) phi(mathbf{y})] mathrm{d} S(mathbf{y}), quad forall mathbf{x} in V
\
end {equation}
\ where G(x,y) is the Green’s function given by \
begin {equation}
G(mathbf{x}, mathbf{y})=frac{1}{2 pi} ln left(frac{1}{r}right)
 end {equation}
\ for 2D problems, and \
begin {equation}
F(mathbf{x}, mathbf{y})=frac{partial G(mathbf{x}, mathbf{y})}{partial n(mathbf{y})}=frac{1}{2 pi r} frac{partial r}{partial n}
 end {equation}
with r being the distance between the collocation point x and field point y (Fig.1).
Letting x/S (Fig. 1), we obtain the following boundary integral equation [7–9] \
begin {equation}
C(mathbf{x}) phi(mathbf{x})=int_{S}[G(mathbf{x}, mathbf{y}) q(mathbf{y})-F(mathbf{x}, mathbf{y}) phi(mathbf{y})] mathrm{d} S(mathbf{y}), quad forall mathbf{x} in S
 end {equation}
in which the second-term on the right-hand side is a singular integral of the Cauchy-principal value (CPV) type, and the coefficient \
begin {equation}
C(mathbf{x})=-int_{S} F(mathbf{x}, mathbf{y}) mathrm{d} S(mathbf{y})
 end {equation}
which is also a CPV integral. If the boundary S is smooth at the collocation point x, we simply have C(x)Z1/2. Note that by substituting (8) into (7), we can obtain the following weakly-singular form of the BIE (see, e.g. Refs. [22,23]) for potential problems \
begin {equation}
int_{S} F(mathbf{x}, mathbf{y})[phi(mathbf{y})-phi(mathbf{x})] mathrm{d} S(mathbf{y})=int_{S} G(mathbf{x}, mathbf{y}) q(mathbf{y}) mathrm{d} S(mathbf{y}), quad forall mathbf{x} in S
end {equation}
in which no singular integrals exist. One can choose either BIE (7) or (9) in the discretization using the BEM, depending on which form is more convenient.
As an example of discretization schemes, we use constant boundary elements, that is, dividing the boundary S into N line segments (elements) and placing one node on each element (Fig. 2). We obtain the following discretized equation of BIE (7) for node i [7–9] \
begin {equation}
frac{1}{2} phi_{i}=sum_{j=1}^{N}left[g_{i j} q_{j}-f_{i j} phi_{j}right], quad text { for } i=1,2, ldots, N
end {equation}

begin{wrapfigure}{l}{linewidth}
vspace{-20pt}
begin{center}
includegraphics[width=linewidth]{./fig_2.png}
end{center}
vspace{-20pt}
caption{Figure}
end{wrapfigure}

indent \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \
where fj and qj (jZ1,2,.,N) are the nodal values off and q on element DSj (Fig. 2), respectively, and the coefficients are given by
begin {equation}
phi(mathbf{x})=overline{phi}(mathbf{x}), quad forall mathbf{x} in S_{1}
 end {equation}
with the collocation point x being placed at node i. In matrix form, Eq. (10) is written as \
%               HERE IS PROBLEM!    \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\
In the conventional BEM approach, a standard linear system of equations is formed as follows by applying the boundary condition (either Eq. (2) or (3)) at each node and switching the columns in the two matrices in Eq.(12) \
begin {equation}
begin{array}{l}{left[ begin{array}{cccc}{a_{11}} & {a_{12}} & {cdots} & {a_{1 N}} \ {a_{21}} & {a_{22}} & {cdots} & {a_{2 N}} \ {vdots} & {vdots} & {ddots} & {vdots} \ {a_{N 1}} & {a_{N 2}} & {cdots} & {a_{N N}}end{array}right] left{begin{array}{c}{lambda_{1}} \ {lambda_{2}} \ {vdots} \ {lambda_{N}}end{array}right}=left{begin{array}{c}{b_{1}} \ {b_{2}} \ {vdots} \ {b_{N}}end{array}right},} \ \ {mathbf{A} lambda=mathbf{b}}end{array}
 end {equation}
where A is the coefficient matrix, l the unknown vector and b the known right-hand side vector. Obviously, the construction of matrix A requires O(N2) operations using the two expressions in Eq. (11) and the size of the required memory for storing A is also O(N2) since A is in general a non-symmetric and dense matrix. The solution of system in Eq. (13) using direct solvers such as Gauss elimination is even worse, requiring O(N3) operations because of this general matrix. That is why the conventional BEM approach for solving the BIE’s is in general slow and inefficient for large-scale problems, despite its robustness in the meshing stage as compared with other domain based methods.

end{multicols}


end{document}