A Fast Multipole Accelerated Singular Boundary Method for Potential Problems

W.Chen,C.J.Liuand Y.Gu

A Fast Multipole Accelerated Singular Boundary Method for Potential Problems

W.Chen1,2,C.J.Liu1and Y.Gu2,3

The singular boundary method(SBM)is a recently-developed meshless boundary collocation method.This method overcomes the well-known f i ctitious boundary issue associated with the method of fundamental solutions(MFS)while remaining the merits of the later of being truly meshless,integral-free,and easy-to-program.Similar to the MFS,this method,however,produces dense and unsymmetrical coeff i cient matrix,which although much smaller in size compared withdomaindiscretizationmethods,requiresO(N2)operationsintheiterativesolution of the resulting algebraic system of equations.To remedy this bottleneck problem for its application to large-scale problems,this paper makes the f i rst attempt to develop a fast multipole SBM(FM-SBM)formulation for two-dimensional(2D)potential problems.The proposed strategy can solve large-scale problems with several millions boundary discretization nodes on a desktop computer.

Singular boundary method,method of fundamental solutions,fast multipole method,meshless boundary collocation method,potential problem.

1 Introduction

The singular boundary method(SBM)[Chen and Wang(2010)]is a recentlydeveloped meshless boundary collocation method[Atluri et al.(2006);Dong and Atluri(2011),(2012);Zhang et al.(2002);Zhang et al.(2003)]and has since been successfully applied to a variety of physical problems,such as potential[Chen and Gu(2012);Gu and Chen(2013);Gu et al.(2012a)]and elastostatic problems[Gu et al.(2014);Gu et al.(2012b);Gu et al.(2011)].This method inherits the merits of the method of fundamental solutions(MFS)[Chen(1995);Fairweather andKarageorghis(1998);Sarler(2009);Sarler and Liu(2013);Young et al.(2007)]and is truly free of mesh and integration.The SBM avoids the perplexing f i ctitious boundary issue associated with the traditional MFS.However,the application of this method has so far been limited to model relatively small scale problems.This is mainly because the SBM produces dense interpolation matrix that,although smaller in size compared with the domain discretization methods such as the FEM,requires O(N2)operations in the iterative solution of its resulting algebraic equations.The SBM solution of large-scale problems is computationally expensive.

Inspired by the successful application of the fast multipole method(FMM)[Greengard and Rokhlin(1987);Rokhlin(1985)]in signif i cantly reducing computing costs of the desnse matrix equations from boundary integraltion methods(BIEs)[Liu and Nishimura(2006);Liu et al.(2005);Nishimura(2002);Wang and Yao(2004)],this study makes the f i rst attempt to apply the FMM to accelerating the solutions of the SBM discretization equations.Applying the FMM to the SBM can potentially reduce both the operation count and memory requirement to O(N).This will be a signif i cant improvement to the SBM in extending its applications.

This paper demonstrates that the proposed FM-SBM strategy is capable of solving large-scale problems.The FM-SBM formulations and its numerical implementation will be discussed in details.Numerical examples with up to 4,000,000 collocation points are successfully solved on a Core II PC using the developed FM-SBM code.

A brief outline of the rest of this paper is as follows:Section 2 reviews the conventional SBM formulations for potential problems.Section 3 provides detailed formulations and algorithms used in the FM-SBM,followed in Section 4 two benchmark numerical examples are well studied to illustrate the accuracy and eff i ciency of the proposed strategy.Finally,the conclusions and remarks are provided in Section 5.

2 Formulations of the conventional SBM

Without loss of generality,we consider the following Laplace equation governing potential problems in a 2D domain Ω(see Fig.1)of the coincident source-collocation point for mathematically well-posed problems.Our numerical experiments illustrate that the origin intensity factors do exist and have f i nite values,only depending on the distribution of boundary nodes and their respective boundary conditions.We refer readers to Refs.[Gu et al.(2014);Gu et al.(2011)]for further detailed derivation of the aforementioned origin intensity factors.

3 The fast multipole method for SBM

The key idea of the FMM is to convert the point-to-point interactions to cell-to-cell interactions,in which the procedure is eff i ciently carried out with tree structures.An iterative solver,such as GMRES,is employed to reduce the computational complexity of the matrix-vector multiplication from O(N2)to O(N).Using the FMM,both the solution time and the memory requirements can be reduced to order of O(N)[Rokhlin(1985)].

This section will present a FMM formulation for the SBM,which is similar to that for the BEM,except for the calculation of the kernel expansions where integrals are replaced by summations.For completeness,all the required formulas are provided.The detailed derivation of some of these formulas can be found in Refs.[Liu and Nishimura(2006);Liu et al.(2005)].

For convenience,complex notations are used to represent the collocation point x and the source point s by z0=x0+iy0and z=x1+iy1,respectively(see Fig.2).The fundamental solution(6)can then be rewritten as

in which the real part gives the values of fundamental solution.And the kernels A and B shown blow related to the G kernel,can be rewritten as

with n=n1+in2,in complex notation.

3.1 Multipole moments

Suppose that zcis a point close to source point z,that is,inequality|z−zc|≤is valid(see Fig.2).Using the subtracting and adding-back technique,

and only need to be computed once.

In addition,for the sum on the right-hand side of Eq.(5),we can also derive the following multipole expansion

by using the same moments given in Eq.(14).

In order to calculate the diagonal coeff i cients of Neumann boundary Eq.(5),the sum∑jljB(z0,zj)/lishould be calculated.By using Eqs.(10)and(12),we can obtain the following multipole equation

are another group of moments about zc.

3.2 M2M translation

Suppose that the point zcis moved to a new location zc0(see Fig.2),we obtain the following equations

When applying the binomial formula and rearranging the terms to Eqs.(18)and(19),respectively,we can further get

Note that Eqs.(20)and(21)translate the moments centered at zcto the moments centered at zc0,which are called M2M translation.

3.3 Local expansions and M2L translation

SupposezLisapointclosetothecollocationpointz0,thatis,theinequality|z0−zL|＜＜is valid(see Fig.2).Using the multipole expansions Eqs.(13)and(16),we can get

respectively.

Expanding f(z0)and g(z0)about zLusing the Taylor series expansion,we can derive the following local expansions

where the coeff i cients are given by the following equations

which are called M2L translation.

Using Eq.(26),the local expansion for the A kernel sum can be naturally obtained by

3.4 L2L translation

Suppose the point for local expansion is moved from zLto zL0(see Fig.2),using

Eq.(24),we obtain the following local equation

where the new coeff i cients are given by the following L2L translation

Note that Eqs.(24)and(25)have the same form,and the L2L translations for B kernel sum can be easily obtained by replacing Ll(zL0)and Lm(zL)with(zL0)andrespectively.

It is noted that the M2M and L2L translations for A kernel sum and G kernel sum remain the same.

3.5 FMM formulation for the SBM

The algorithms or procedures for the FM-SBM are similar to those for the FMBEM[Liu et al.(2005);Nishimura(2002)].The main steps in the FM-SBM can be described as follows:

Step 1.Discretization.Discretize the boundary S in the same manner as in the conventional SBM,that is,place a certain number of collocation points,which are also the source points at the same time,on S as shown in Fig.1.associated with cells in the interaction list,and the latter by using L2L translations which shifts the local expansions from the centroid of C,s parent to the centroid of cellC(see Fig.2).We repeat these procedures starting from level 2 and tracing the tree structure downward to all leaves.For a cell at level 2,the coef fi cients of the local expansions are only involved M2L translation.

Step 5.Evaluation of the sum in Eq.(4)or(5).LetC be a leaf to which the collocation point z0belongs,and the values of αjbe given.We calculate the contributions from source points in leafC and its adjacent cells directly as in conventional SBM.We also compute contributions from all other cells by using the coef fi cients of the local expansions associated withC which have been computed in Step 4.

Step 6.First,we let αjequals 1 to compute the diagonal coef fi cients of the coeffi cient matrix in Neumann boundary equation(5)using Step 5.Then,we update the unknown vector in the system produced by using Eqs.(4)or(5),and go back to Step 3 for the matrix and unknown vector multiplication until the solution converges within a given residue using generalized minimal residual(GMRES)solver.

3.6 Preconditioners for the FM-SBM

Similar to FM-BEM,a good pre-conditioner for FM-SBM is crucial for its convergence and computing eff i ciency.Since the conditioning number of the SBM is almost the same as that of the BEM,the present FM-SBM directly utilizes the block diagonal pre-conditioner proposed for FM-BEM[Wang and Yao(2004)].The preconditioner is based on leaves of the hierarchical tree,which is formed on each leaf using direct evaluation of the fundamental solution with the collocation points within that leaf and their corresponding source points.This pre-conditioner will be ineff i cient if the number of collocation points in a leaf is large.

4 Numerical examples

To verify the feasibility and eff i ciency of the proposed FM-SBM,two numerical examples are examined in which the results are compared with the corresponding exact solutions and those obtained using conventional SBM.To assess the accuracy of the FM-SBM,the numerical results will be compared with exact solutions in terms of the relative error def i ned by

Table 1:Relative errors of computed potential and normal derivative.

The comparisons of CPU time and memory requirement between FM-SBM and conventional SBM are illustrated in Figs.5 and 6,respectively.To show the order of complexity of FM-SBM more clearly,logarithmic scale was used in Figs.5 and 6.The vertical axis of Fig.5 denotes the 10-base logarithmic of CPU time,simplifi ed as log(t),and the horizontal axis represents the 10-base logarithmic of number of collocation points,simpli fi ed as log(N).Similarly,the vertical axis of Fig.6 denotes the 10-base logarithmic of memory requirement,simpli fi ed as log(m),and the horizontal axis represents the 10-base logarithmic of number of collocation points,simpli fi ed as log(N).These fi gures show a signi fi cant edge of the present FM-SBM in the savings of CPU time and memory,which are observed both approximately proportional to the number of collocation points.The results verify that both multiplication operations and storage costs of FM-SBM are O(N).This shows that the FM-SBM has great promise to handle large-scale problems.

Fig.7 show the iteration number and maximum number of boundary points in a leaf node vary with the number of collocation points,respectively,as the residue for iteration and the maximum of points allowed in a leaf of the tree structure were fi xed.It can be concluded that,in general,the iteration number and the maximum number of boundary increase with the number of collocation points.

Figure 5:Comparison of the CPU time for the problem of the amoeba-like domain with ten holes.

Figure 6:Comparison of the memory requirement for the problem of the amoebalike domain with ten holes.

Figure 7:Relationships of iteration number and maximum number of points in leaves with the number of collocation points.

Figure 8:Iteration number versus residue of iteration under different scale of problems.

The relationship between the iteration number and residue of iteration is given in Fig.8.From this f i gure,we can obtain that,when the residue of iteration gets less,FM-SBM needs more iteration steps to get the resolution.And,this tendency doesn’t matter with the scale of the problem.We can also get that,for the same residue of iteration,when the scale becomes larger,FM-SBM will need more iteration steps.

Table 2:Relative errors of computed potential and normal derivative.

In Table 2,Nleafand Nleveldenote the number of collocation points in a leaf and the number of levels of the tree structure,respectively.Table 2 shows that,when the number of collocation points increases,the number of levels decreases.As we also know,when the number of collocation points increases,it will need more CPU time to calculate the interaction of nodes in a leaf.So,as illustrated in Table 2,there will be a better value of the number of collocation points which is 400 in this example to get less CPU time.

4.2 A square plate with numerous circular holes

Here,the thermal conductivity of a square plate with uniformly distributed circular holes,or a perforated plate,is studied using the FM-SBM(see Fig.9).Several models of the plate,with increasing number of holes,are considered.Each model of the plate contains a total of 2m×2m holes,with m=2,3,...,28,30.The radii of the holes are the same and equal to 4.0/(2m+2).Each hole is discretized with 100 collocation points and the outer boundary of the plate is discretized with 400mcollocation points.

In this example,the analytical solution is given by

The Dirichlet boundary conditions are specif i ed on the outer boundary,and the Numann boundary condition are specif i ed on the remaining boundaries.For FMSBM,the number of terms equals 20 for multipole and local expansions,and the maximum of collocation points allowed in a leaf is set to 400.The residue for iteration is prescribed to 1×10−6,and the initial guess is zero.

Figure 9:A model of the square plate with 400 uniformly distributed circular holes

Table 3:Relative error of computed potential and normal derivative for the square plate with circular holes.

Table 3 illustrates the relative errors of u and q by both the FM-SBM and conventional SBM.As shown in this table,the conventional SBM can only provide accurate numerical solution when the number of circular holes less than 100,after which the method becomes infeasible due to memory limitations.In sharp contrast,the proposed FM-SBM is still workable when the number of circular holes is as large as 3,600,with the largest relative error less than 3E-5.

Figure 10:Comparison of the CPU time for the square plate with circular holes.

Figure 11:Comparison of the memory for the square plate with circular holes.

The relationships of iteration number and maximum number of points in leaves with the number of collocation points,is shown in Fig.12.From this f i gure,we can conclude that the iteration number and the maximum number of boundary increase with the number of collocation points,generally.

Figure 12:Relationships of iteration number and maximum number of points in leaves with the number of collocation points.

5 Conclusions

In this paper,a fast multipole accelerated singular boundary method is presented for 2D potential problems.Complete formulations and implementation details for FM-SBM are provided.Two numerical examples demonstrate the eff i ciencies of the present strategy for solving large-scale problems.

The convergence issue in the FM-SBM requires further investigations,and a better pre-conditioner is required to devise.Furthermore,the algorithms of the FM-SBM can also be improved.In this work,the calculation of the potential and its derivatives were done by direct calculation using SBM interpolation formula,which may further be speeded up by applying FMM,when the number of data points is comparable to that of collocation points.All these issues merit further investigations based on this study.

Acknowledgement: The work described in this paper was supported by the National Basic Research Program of China(973 Project No.2010CB832702),the National Science Funds for Distinguished Young Scholars of China(No.11125208),andtheNationalScienceFundsofChina(No.11402075,11302069),andtheChina Postdoctoral Science Foundation(No.2015M570572).

Atluri,S.N.;Liu,H.T.;Han,Z.D.(2006):Meshless local Petrov-Galerkin(MLPG)mixed collocation method for elasticity problems.CMES:Computer Modeling in Engineering&Sciences,vol.4,no.3,pp.141-152.

Chen,C.S.(1995):The method of fundamental solutions for non-linear thermal explosions.Commun.Numer.Methods.Eng.,vol.11,no.8,pp.675-681.

Chen,W.;Gu,Y.(2012):An improved formulation of singular boundary method.Adv.Appl.Math.Mech.,vol.4,no.5,pp.543-558.

Chen,W.;Wang,F.Z.(2010):A method of fundamental solutions without f i ctitious boundary.Eng.Anal.Bound.Elem.,vol.34,no.5,pp.530-532.

Dong,L.;Atluri,S.N.(2011):Development of T-Trefftz four-node quadrilateral and Voronoi Cell Finite Elements for macro-µmechanical modeling of solids.CMES:Computer Modeling in Engineering&Sciences,vol.81,no.1,pp.69-118.

Dong,L.;Atluri,S.N.(2012):A simple multi-source-point Trefftz method for solvingdirect/InverseSHMproblemsofplaneelasticityinarbitrarymultiply-connecte domains.CMES:Computer Modeling in Engineering&Sciences,vol.85,no.1,pp.1-44.

Fairweather,G.;Karageorghis,A.(1998):The method of fundamental solutions for elliptic boundary value problems.Adv.Comput.Math.,vol.9,no.1,pp.69-95.

Greengard,L.;Rokhlin,V.(1987):A fast algorithm for particle simulations.J.Comput.Phys.,vol.73,no.2,pp.325-348.

Gu,Y.;Chen,W.(2013):Inf i nite domain potential problems by a new formulation of singular boundary method.Appl.Math.Model.,vol.37,no.4,pp.1638-1651.

Gu,Y.;Chen,W.;He,X.-Q.(2012a):Singular boundary method for steady-state heat conduction in three dimensional general anisotropic media.Int.J.Heat.Mass.Transf.,vol.55,no.17–18,pp.4837-4848.

Gu,Y.;Chen,W.;He,X.(2014):Improved singular boundary method for elasticity problems.Comput.Struct.,vol.135,no.0,pp.73-82.

Gu,Y.;Chen,W.;He,X.Q.(2012b):Domain-Decomposition Singular Boundary Method for Stress Analysis in Multi-Layered Elastic Materials.CMC:Computers Materials&Continua,vol.29,no.2,pp.129-154.

Gu,Y.;Chen,W.;Zhang,C.-Z.(2011):Singular boundary method for solving plane strain elastostatic problems.Int.J.Solids.Struct.,vol.48,no.18,pp.2549-2556.

Liu,Y.J.;Nishimura,N.(2006):The fast multipole boundary element method for potential problems:A tutorial.Engineering Analysis with Boundary Elements,vol.30,no.5,pp.371-381.

Liu,Y.J.;Nishimura,N.;Yao,Z.H.(2005):A fast multipole accelerated method of fundamental solutions for potential problems.Engineering Analysis with Boundary Elements,vol.29,no.11,pp.1016-1024.

Nishimura,N.(2002):Fastmultipoleacceleratedboundaryintegralequationmethods.Applied Mechanics Reviews,vol.55,no.4,pp.299-324.

Rokhlin,V.(1985):Rapid solution of integral equations of classical potential theory.Journal of Computational Physics,vol.60,no.2,pp.187-207.

Sarler,B.(2009):Solution of potential f l ow problems by the modif i ed method of fundamental solutions:Formulations with the single layer and the double layer fundamental solutions.Eng.Anal.Bound.Elem.,vol.33,no.12,pp.1374-1382.

Sarler,B.;Liu,Q.G.(2013):Non-singular method of fundamental solutions for two-dimensional isotropic elasticity problems.CMES:Computer Modeling in Engineering&Sciences,vol.91,no.4,pp.235-266.

Wang,H.T.;Yao,Z.H.(2004):Application of a new fast multipole BEM for simulation of 2D elastic solid with large number of inclusions.Acta Mechanica Sinica,vol.20,no.6,pp.613-622.

Young,D.L.;Chen,K.H.;Chen,J.T.;Kao,J.H.(2007):A modif i ed method of fundamental solutions with source on the boundary for solving Laplace equations with circular and arbitrary domains.CMES:Computer Modeling in Engineering&Sciences,vol.19,no.3,pp.197-221.

Zhang,J.M.;Yao,Z.H.;Li,H.(2002):A hybrid boundary node method.Int.J.Numer.Methods.Eng.,vol.53,no.4,pp.751-763.

Zhang,J.M.;Yao,Z.H.;Tanaka,M.(2003):The meshless regular hybrid boundary node method for 2D linear elasticity.Eng.Anal.Bound.Elem.,vol.27,no.3,pp.259-268.

1State Key Laboratory of Hydrology-Water Resources and Hydraulic Engineering,International Center for Simulation Software in Engineering and Sciences,College of Mechanics and Materials,Hohai University,Nanjing 210098,China.

2Corresponding author.E-mail:chenwen@hhu.edu.cn;guyan1913@163.com

3College of Mathematics,Qingdao University,Qingdao 266071,PR China.