# Mechanical Quadrature Method and Splitting Extrapolation for Solving Dirichlet Boundary Integral Equation of Helmholtz Equation on Polygons

Authors: Hu Li; Yanying Ma
Journal: Journal of Applied Mathematics
 (2014)

---

## Abstract

We study the numerical solution of Helmholtz equation with Dirichlet boundary condition. Based on the potential theory, the problem can be converted into a boundary integral equation. We propose the mechanical quadrature method (MQM) using specific quadrature rule to deal with weakly singular integrals. Denote byhm the mesh width of a curved edge Γm(m=1,…,d) of polygons. Then, the multivariate asymptotic error expansion of MQM accompanied with O(hm3) for all mesh widths hm is obtained. Hence, once discrete equations with coarse meshes are solved in parallel, the higher accuracy order of numerical approximations can be at least O(hmax⁡5) by splitting extrapolation algorithm (SEA). A numerical example is provided to support our theoretical analysis.


---

## Body

## 1. Introduction

We consider Helmholtz equation with Dirichlet boundary condition:(1)Δu(x)+α2u(x)=0,x∈Ω,um(x)=gm(x),x∈Γm(Γ=⋃m=1d),
where Ω⊂R2 is a polygonal domain with the curved boundary Γ=∪m=1dΓm, d>1, and the angels between Γm and Γm+1 are in (0,2π], and the function gm(x)=g(x)|Γm is known on Γm.By the potential theory, the solutions of (1) can be represented as a single-layer potential
(2)u(y)=∫Γ‍Φ(y,x)v(x)dsx,y=(y1,y2)∈Ω,
where x=(x1,x2) and Φ(y,x) is the foundation solution of Helmholtz equation, which is given by
(3)Φ(y,x)=i4H0(1)(α|x-y|).
And H0(1)=J0+iN0 is the Hankel functions of order zero and of the first kind, where
(4)J0(z)=∑n=0∞‍(-1)n(n!)2(z2)2n
for the Bessel function of order zero and
(5)N0(z)=2π(ln⁡z2+γ)J0(z)+2π∑n=1∞‍(∑m=1n‍1m)(-1)n+1(n!)2(z2)2n
for the Neumann function of order zero. And γ=0.57721… is Euler constant.In what follows, in order to analyze properties of the kernel, we decompose the kernel(6)Φ(x,y)=k1(x,y)+k2(x,y),
where k1(x,y)=-1/2πln⁡|x-y| is logarithmic singular function and k2(x,y)=i/4-1/2π(ln⁡(α/2)+γ)+O(|x-y|ln⁡|x-y|) is a smooth function. v(x) is the solution of the following equation:
(7)g(y)=∫Γ‍Φ(y,x)v(x)dsx,y=(y1,y2)∈Γ.
Equation (7) is weakly singular BIE system of the first kind, whose solution exists and is unique as CT≠1 [1], where CT is the logarithmic capacity (i.e., the transfinite diameter). Once v(x) is solved from (7), the function u(y)(y∈Ω) can be calculated by (2).Galerkin and collocation methods [2, 3] are used to solve (7). However, the discrete matrix is full and each element has to calculate the weakly singular integral for collocation methods or the double weakly singular integral for Galerkin methods, which implies that the work calculating discrete matrix is so large as greatly to exceed to solve the discrete equations. When the numerical methods are applied, the accuracy of numerical solutions is lower at singular points [4] and the corresponding numerical results become unreliable, because the condition numbers are very large.In the paper, MQM is proposed to calculate weakly singular integrals by Sidi quadrature rules [5], which makes the calculation of the discrete matrix becomes very simple and straightforward without any singular integrals. The convergence theory of approximations is given by estimating eigenvalues of the discrete matrix and using Anselone’s collectively compact convergent theory [6], which shows that the method retains the optimal convergence order O(hmax⁡3) and possesses the optimal condition number O(hmin⁡-1). Since MQM possesses the multivariate asymptotic expansion of errors, we can construct SEA to obtain the convergence order O(hmax⁡5). Once discrete equations on some coarse meshes are solved in parallel, the accuracy of numerical solutions can be greatly improved by SEA.This paper is organized as follows. Section2 includes the singularity analysis of the integral kernels and the solution. In Section 3, the MQM is described. In Section 4, we can obtain multiparameter asymptotic expansion of errors and SEA is described. In Section 5, a numerical example is provided to verify the theoretical results.


## 2. The Integral Kernels and the Solution of Singularity Analysis

Define boundary integral operators onΓm:
(8)(Kqmvm)(y)=-12π∫Γm‍Φ(x,y)vm(x)dsx,ddddy∈Γq,m,q=1,…,d.
Then (7) can be converted into a matrix operator equation
(9)KV=G,
where K=[Kqm]q,m=1d, V=(v1(x),…,vd(x))T, and G=(g1(y),…,gd(y))T.Assume thatΓm can be described by the parameter mapping:
(10)xm(s)=(x1m(s),x2m(s)):[0,1]⟶Γm,0≤s≤1,
with |xm′(s)|=[(x1m′(s))2+(x2m′(s))2]1/2>0, m=1,…,d. Using Sidi periodic transformation [7]
(11)ψp(τ)=ϑp(τ)ϑp(1):[0,1]⟶[0,1],p∈N,
where ϑp(τ)=∫0τ‍(sin(πρ))pdρ, then operator (8) is also converted into integral operator on [0,1]:
(12)(B¯q0v¯q)(t)=-12π∫01‍ln⁡|2e-1/2sin(π(t-τ))|v¯q(τ)dτ,(B¯0qv¯q)(t)=-12π∫01‍(ln⁡|x¯q(t)-x¯q(τ)2e-1/2sin(π(t-τ))|-πi2dddddddddd+ln⁡α2+γ|x¯q(t)-x¯q(τ)2e-1/2sin(π(t-τ))|)v¯q(τ)dτ,(K¯qmv¯m)(t)=∫01‍Φ(x¯q(t),x¯m(τ))v¯m(τ)dτ,q≠m,
where x¯q(t)=xq(ψp(t)) and v¯m(τ)=vm(ψp(τ))|xm′(ψp(τ))|ψp′. Then (9) can be rewritten as
(13)(B¯0+B¯1+K¯)V¯=G¯,
where B¯0=diag⁡(B¯10,…,B¯d0), B¯1=diag⁡(B¯01,…,B¯0d), K¯=[K¯qm]q,m=1d, V¯=(v¯1(τ),…,v¯d(τ))T, G¯=(G¯1(t),…,G¯d(t))T, and G¯m(t)=gm(ψp(t)).Letbq0(t,τ), b0q(t,τ), and kqm(t,τ) be the kernels of the integral operators B¯q0, B¯0q, and K¯qm, respectively. Then, the following results hold.(1)
bq0(t,τ) is a logarithmic singular function on t∈[0,1] and τ∈[0,1].(2)
b0q(t,τ) is a continuous function on t∈[0,1] and τ∈[0,1].(3)
ForΓq∩Γm=∅(|q-m|≠1ord-1), kqm(t,τ) is a continuous function on t∈[0,1] and τ∈[0,1].(4)
The solutionv¯q(t) is a smooth function under (11).Lemma 1.
IfΓq∩Γm≠∅(|q-m|=1 or d-1) and p≥2, then k¯qm(t,τ)(=sinp(πt)kqm(t,τ)) and (∂n/∂tn)k¯(t,τ)(n=1,2) are continuous on [0,1]2.Proof.
Without loss of generality, we assume that the origin(0,0)=Γm∩Γq of coordinates is a vertex with an interior angle θm+1. By (6), kqm can be expressed by
(14)kqm(t,τ)=-12πln⁡[a02(t)+a1(τ)2-2a0(t)a1(τ)cos⁡θm+1]+k2(xq,xm),
where a0(t)=|xq(ψp(t))| and a1(τ)=|xm(ψp(τ))|. We assume that a0(0)=a1(0)=0; namely, kqm has the logarithmic singularity at points (0,0) and is continuously differentiable in [0,1]2∖(0,0). Consider
(15)kqm(t,τ)=-12πln⁡(a02(t)+a1(τ)2)-12π×ln⁡[1-2a0(t)a1(τ)cos⁡θm+1a02(t)+a1(τ)2]+k2(xq,xm).
We have
(16)|2a0(t)a1(τ)cos⁡θm+1a02(t)+a1(τ)2|≤|cos⁡θm+1|<1.
So the singularity of kqm(t,τ) occurs in (1/2π)ln⁡(a02(t)+a1(τ)2). In order to remove the singularity, we define
(17)ς(t,τ)=sinp(πτ)ψp(1)ln⁡(a02(t)+a1(τ)2).
Since ψp′ have p order zero, we have the error estimate
(18)|ς(t,τ)|=o(ɛp)|ln⁡ɛ|⟶0,asɛ⟶0,
which shows that ς(t,τ) is bounded and continuous on (C2[0,1])2 when ɛ/2≤t, τ≤ɛ for any ɛ>0.
Similarly, we can prove(∂/∂t)ς(t,τ) and (∂2/∂t2)ς(t,τ) are continuous on (C[0,1])2, which can be obtained by
(19)|∂∂tς(t,τ)|≤|sinp(πτ)ψp(1)2a0(t)|xm+1′(ψp(t))||ψp′(t)|a02(t)+a12(τ)|=O(ɛp),∀t,τ∈[ɛ2,ɛ],
and (∂2/∂t2)ς(t,τ)=O(ɛp-1). As above, (∂n/∂tn)ς(t,τ)(n=0,1,…,p) are continuous on (C[0,1])2. We can prove that (∂n/∂tn)k¯qm(t,τ)(n=0,1,…,p) are also continuous on (C[0,1])2. The proof of Lemma 1 is completed.


## 3. Mechanical Quadrature Method

Lethm=1/nm(nm∈N,m=1,…,d) be mesh widths and let tmj=(j-1/2)hm(j=1,…,nm) be nodes. Since an integral operator W is continuous, by the trapezoidal or the midpoint rule [8], we can construct Nyström approximate operator Wh for the integral operator W, defined by
(20)(Whv¯m)(t)=hm∑j=1nm‍w(t,tmj)v¯m(tmj),
and the error
(21)(Whv¯m)(t)-(Wv¯m)(t)=O(hm2l),l∈N.
Since B¯q0 has the singularities on [0,1], by Sidi quadrature formula [5], we get the following approximations B¯q0h of integral operator B¯q0:
(22)(B¯q0hv¯m)(t)=-12πhm∑j=1,t≠tmjnm‍ln⁡|2e-1/2sinπ(t-tmj)|×v¯m(tmj)-hm2πln⁡|2πe-1/2hm(2π)|v¯m(t),
and the error
(23)(B¯q0hv¯m)(t)-(B¯q0v¯m)(t)=-2π∑μ=12l-1‍ζ′(-2μ)(2μ)![v¯m(t)]2μhm2μ+1+O(hm2l),
where ζ′(t) is the derivative of the Riemann zeta function.Then, we can get the approximate equation of (13):
(24)(B¯0h+B¯1h+K¯h)V¯h=G¯h,
where V¯h=(V¯1h,…,V¯dh), V¯mh=(v¯m(t1),…,v¯m(tnm)), B¯0h=diag⁡(B¯10h,…,B¯d0h), B¯q0h=[bq0(tj,τi)]j,i=1nq, B¯1h=diag⁡(B¯01h,…,B¯0dh), B¯0q=[b0q(tj,τi)]j,inq, K¯h=[K¯qmh]q,m=1d, K¯qmh=[kqm(tj,τi)]j,i=1nq,nm, G¯h=(G¯1h,…,G¯dh)T, and G¯qh=(gq(t1),…,gq(tnq)). Obviously, (24) is a system of linear equations with n(=∑m=1d‍nm) unknowns. Once V¯h is solved by (24), the solution u(y)(y∈Ω) can be computed by
(25)uh(y)=∑m=1d‍hm∑j=1nm‍Φ(y,x(tmj))v¯m(tmj).B¯q0h is symmetric circular matrix and has the form of
(26)B¯q0h=circulate(-hm2πln⁡|2e-1/2hmπ|,…,hhhhhhihhhh-hm2πln⁡|2e-1/2sin(π(nm-1)hm)|).Lemma 2 (see [9]).
(1) There exists a positivec1>0 so that the eigenvalues λβ(β=1,…,nm) of B¯q0h satisfy c1>λβ>1/(2πnm). (2) The condition number of B¯q0h is O(nm). (3) B¯q0 is invertible, and (B¯q0h)-1 is uniformly bounded with the spectral norm ∥(B¯q0h)-1∥=O(nm).Based on Lemma2, we immediately get the following corollary.Corollary 3.
(1)B¯0h is invertible, and (B¯0h)-1 is uniformly bounded with the spectral norm ∥(B¯0h)-1∥=O(n0). (2) The condition number of B¯0h is O(n0), where n0=max⁡m=1dnm.From Corollary3, we know that (24) is equivalent to
(27)(Eh+(B¯0h)-1(B¯1h+K¯h))V¯h=(B¯0h)-1Gh,
where Eh denotes the unit matrix.Now we give the following definitions to discuss the approximate convergence in (27).Define a subspaceC0[0,1]={v(t)∈C[0,1]:v(t)/sin2(πt)∈C[0,1]} of the space C[0,1] with a norm ∥v∥*=max⁡0≤t≤1|v(t)/sin2(πt)|. Let Smh⊂C0[0,1] be a piecewise linear function subspace with base points {ti}i=0nm-1, and let ej(t)j=0,…,nm-1 be base function satisfying ej(ti)=δji. Define a prolongation operator Phm:Rhm→Shm satisfying
(28)Phmω=∑j=0nm-1‍ωmej(t),∀ω=(ωm1,…,ωmnm)∈Rnm.
Define a restricted operator Rhm:C0[0,1]→Rnm satisfying
(29)Rhmvm=(vm(tm1),…,vm(tmnm))∈Rnm,kkkkjkkkkkkkkkkk∀vm∈C0[0,1].Replacing(B^0h)-1(B^1h+K^h)=Ph(B¯0h)-1Rh(B¯1h+K¯h), V^h=PhV¯h, and G^=PhG¯h, we construct an operator equation
(30)(Eh+(B^0h)-1(B^1h+K^h))V^h=(B^0h)-1G^h,
where Ph=diag⁡(P1,…,Pd) and Rh=diag⁡(R1,…,Rd). Obviously, if V^h is the solution of (30), then RhV¯h must be the solution of (27); conversely, if V¯h is the solution of (27), then PhV¯h must be the solution of (30). In order to prove the convergence of MQM, we give the following lemma.Lemma 4 (see [9]).
The operator sequence{Pq(B¯q0h)-1RqB¯q0:C2[0,1]→C[0,1]} is uniformly bounded and convergent to the embedding operator I.Corollary 5.
Let the Nyström approximationWh be defined by (20). (1) For Γq∩Γm=∅, one has
(31)Pm(B¯m0)-1RmWh→c.c(B¯m0)-1W,inC[0,1]⟶C[0,1].
(2) For Γq∩Γm≠∅, one has
(32)Pm(B¯m0)-1RmK~qmh→c.c(B¯m0)-1K~qm,ddddddddidinC[0,1]⟶C[0,1],
where K~qmh is the Nyström approximation of integral operator K~qm with the kernel k¯qm(t,τ)=sinp(πt)kqm(t,τ) and →c.c denotes the collectively compact convergence.Theorem 6.
LetΓ=∪m=1dΓm with CT≠1 and Γm(m=1,…,d) be smooth curve. The operator sequence {Ph(B¯0h)-1Rh(B¯1h+K¯h)} is collectively convergent to (B¯0)-1(B¯1+K¯) in V=(C0[0,1])d:
(33)Ph(B¯0h)-1Rh(B¯1h+K¯h)→c.c(B¯0)-1(B¯1+K¯).Proof.
LetΘ={υ:∥υ∥≤1,υ∈V} be a unit ball. H={h1,h2,…} is the grid step sequence, where hn={h1(n),h2(n),…,hd(n)} with n→∞, max1≤m≤dhm(n)→0. We arbitrarily take a sequence {Zh}h∈H in space Θ, where Zh=(Zh1,Zh2,…,Zhd)T with ∥Zhm/sin2(πt)∥∞≤1, m=1,…,d. We firstly consider the first component of Ph(B¯0h)-1Rh(B¯1h+K¯h)Zh:
(34)∑m=1d‍P1(B¯10h)-1R1(B¯01h+K¯1mh)Zhm.
For Γm∩Γ1=∅, by Corollary 5, P1(B¯10h)-1R1(B¯01h+K¯1mh)→c.cB¯10-1(B¯01+K¯1m) in C[0,1], and there exists a convergent subsequence in {P1(B¯10h)-1R1(B¯01h+K¯1mh)Zhm}. For Γm∩Γ1≠∅, we have
(35)∥P1(B¯10h)-1R1K¯1mhZhm∥0,0=∥P1(B¯10h)-1R1K~1mh(Zhmsin2(πt))∥0,0≤∥P1(B¯10h)-1R1B¯10∥0,2∥B¯10-1K~1mh∥2,0∥Zhmsin2(πt)∥∞.
By Lemma 4, there exists a convergent subsequence in {P1(B¯10h)-1R1(B¯01h+K¯1mh)Zhm}. Based on the above two cases it is proved that there exists an infinite subsequence H(1)⊂H such that the first component converges. Similarly, it can be concluded that there exists an infinite subsequence H(d)⊂H(d-1)⊂⋯⊂H(1)⊂H such that Ph(B¯0h)-1Rh(B¯1h+K¯h)Zh converges. Therefore, Ph(B¯0h)-1Rh(B¯1h+K¯h) is collectively compact convergent sequence, and Ph(B¯0h)-1Rh(B¯1h+K¯h)→P(B¯0)-1(B¯1+K¯), where →P shows the point convergence. We complete the proof.For the stability of MQM, we have the following corollary.Corollary 7.
LetΓ=∪m=1dΓm with CΓ≠1, let Γm(m=1,…,d) be smooth curve, let B¯0h, B¯1h, and K¯h be the discrete matrices defined by (20) and (22), respectively, and let λi(i=1,…,n) be the eigenvalues of discrete matrix Mh=B¯0h+B¯1h+K¯h. Then, there exists the bound of condition number
(36)Cond(Mh)=max⁡1≤i≤n⁡|λi(Mh)|max⁡1≤i≤n⁡|λi(Mh)|=O(hmin⁡-1),
where hmin⁡=min⁡m=1dhm, hm=1/nm, is the mesh step size of a curved edge Γm.


## 4. Multiparameter Asymptotic Expansion of Errors and SEA

In this section, we derive the multivariate asymptotic expansion of solution errors and describe SEA. We first provide the main result.Theorem 8.
LetΓ=∪m=1dΓm with CΓ≠1. There exists a vector function Φ=(ϕ1,…,ϕd)T independent of h=(h1,…,hd) such that the following multiparameter asymptotic expansion holds at nodes:
(37)V¯h-V¯=diag⁡(h13,…,hdd)Φ+o(h03),h0=max⁡1≤m≤d⁡hm,
where V¯h∈Sd, V¯∈C3[0,1].Proof.
By (20) and (22), there exists the asymptotic expansion
(38)(G¯h-G¯)=diag⁡(h13,…,hdd)PhRhω+o(h03),
where h=max⁡m=1dhm and ω=(ω1,…,ω)T with ωm=-ζ′(-2)(Gm(t))′′/π.
Using (13), (21), (23), and (24), we obtain
(39)(B¯0h+B¯1h+K¯h)Rh(V¯h-V¯)=V¯h-Ph(B¯0h+B¯1h+K¯h)RhV¯=V¯h-[(B¯0+B¯1+K¯)RhV¯ddddid-diag⁡(h13,…,hdd)PhRhγ]+o(h03)=(G¯h-G¯)+diag⁡(h13,…,hdd)PhRhγ+o(h03)=diag⁡(h13,…,hdd)PhRhη+o(h03),
where γ=(γ1,…,γd)T with γm=-ζ′(-2)(x′(t)ψ′(t)v(t))′′/π and η=(η1,…,ηd)T with ηm=ωm+γm. From Corollary 3, we can obtain
(40)(E+Ph(B¯0h)-1Rh(B¯1h+K¯h))(V¯h-V¯)=diag⁡(h13,…,hdd)(B¯0h)-1PhRhη+o(h03).
Constructing the auxiliary equation
(41)(E+(B¯0)-1(B¯1+K¯))Φ=(B¯0)-1η
and its approximate equation
(42)(E+Ph(B¯0h)-1Rh(B¯1h+K¯h))Φh=(B¯0h)-1PhRhη
and substituting (42) into (40), we obtain
(43)(E+Ph(B¯0h)-1Rh(B¯1h+K¯h))×(V¯h-V¯-diag⁡(h13,…,hdd)Φh)=o(h03).
Since (E+Ph(B¯0h)-1Rh(B¯1h+K¯h)) is uniformly bounded from Theorem 6, we get
(44)V¯h-V¯-diag⁡(h13,…,hdd)Φh=o(h03).
Replacing Φh with Φ, we can complete the proof.The multiparameter asymptotic expansion (37) means that SEA can be applied to solve (7); that is, higher order accuracy o(h03) at coarse grid points can be obtained by solving some discrete equations in parallel. The process of SEA is as follows [10].Step 1.
Takeh(0)=(h1(0),…,hd(0)) and h(m)=(h1(0),…,hm(0)/2,…,hd(0)), and solve (24) under mesh parameters h(m) in parallel to get the numerical solutions V¯h(0)(tmj) and V¯h(m)(tmj), m=1,…,d, j=1,…,nm.Step 2.
Computeuh(0)(y) and uh(m)(y)(y∈Ω), by (25), V¯h(0)(tmj), and V¯h(m)(tmj).Step 3.
Compute the extrapolation on the coarse grids as follows:(45)u*(y)=87[∑m=1d‍uh(m)(y)-(d-78)uh(0)(y)].Step 4.
Compute a posteriori error estimate on the coarse grids as follows:(46)|u(y)-1d∑m=1d‍uh(m)(y)|≤|u(y)-[∑m=1d‍uh(m)(y)-(d-78)uh(0)(y)]|+(87d-1)|1d∑m=1d‍uh(m)(y)-uh(0)(y)|≤(87d-1)|1d∑m=1d‍uh(m)(y)-uh(0)(y)|+o(h03).
In the actual calculation process, a posteriori error estimate is immediately used to verify the calculation accuracy.


## 5. Numerical Example

In this section, we carry out a numerical example for the Helmholtz equation by MQM and SEA, in order to verify the error and stability analysis in the previous sections. LetRe(ep) be the real part of eP=|uh(P)-u(P)| and let Im⁡(eP) be the imaginary part of eP=|uh(P)-u(P)|, where P denotes a point. Post-error and SEA-error denote the a posteriori error and the error after SEA once, respectively.Example 1.
Consider Helmholtz equation withα=2 on a plate domain Ω. We describe the boundary Γ=∪m4Γm with Γ1={(x1,x2)∣0≤x1≤1,x2=0}, Γ2={(x1,x2)∣0≤x2≤1,x1=1}, Γ3={(x1,x2)∣0≤x1≤1,x2=1}, and Γ4={(x1,x2)∣0≤x2≤1,x1=0}. Dirichlet boundary conditions corresponding to the analytical solution ei×(x1+x2) were applied to the boundary. We compute the numerical solution uh with P=A(0.5,0.5) and P=B(0.6,0.6) using φ3(t). The numerical results are listed in Tables 1 and 2. From Table 1, we can know that the convergence rates of uh are O(hmax⁡3) for MQM and that the convergence rates of uh are at least O(hmax⁡5) for SEA. From Table 2, we can see Cond|(2k+1,2k+1,2k+1,2k+1)/Cond|(2k,2k,2k,2k)≈2(k=3,4,5) to indicate Corollary 7. It verifies the stability of convergent theory for MQM. Those results coincide with the theoretical analysis made.Table 1
The numerical results.
(
n
1
,
n
2
,
n
3
,
n
4
)
Re
(
e
A
)
Im
⁡
(
e
A
)
Re
(
e
B
)
Im
⁡
(
e
B
)
(8, 8, 8, 8)
3.943
e
-
4
2.459
e
-
4
6.255
e
-
4
2.832
e
-
4
(16, 8, 8, 8)
4.730
e
-
4
2.351
e
-
4
6.756
e
-
4
2.701
e
-
4
(8, 16, 8, 8)
1.400
e
-
4
1.472
e
-
4
3.007
e
-
4
1.688
e
-
4
(8, 8, 16, 8)
1.400
e
-
4
1.472
e
-
4
3.007
e
-
4
1.688
e
-
4
(8, 8, 8, 16)
4.730
e
-
4
2.351
e
-
4
6.756
e
-
4
2.701
e
-
4
Post-error
3.134
e
-
4
1.953
e
-
4
4.902
e
-
4
2.277
e
-
4
SEA-error
6.841
e
-
6
4.125
e
-
6
2.080
e
-
6
8.284
e
-
6
(16, 16, 16, 16)
4.534
e
-
5
2.849
e
-
5
7.934
e
-
5
2.973
e
-
5
(32, 16, 16, 16)
5.438
e
-
5
2.724
e
-
5
8.549
e
-
5
2.788
e
-
5
(16, 32, 16, 16)
1.648
e
-
5
1.729
e
-
5
3.849
e
-
5
1.859
e
-
5
(16, 16, 32, 16)
1.648
e
-
5
1.729
e
-
5
3.849
e
-
5
1.859
e
-
5
(16, 16, 16, 32)
5.438
e
-
5
2.724
e
-
5
8.549
e
-
5
2.788
e
-
5
Post-error
3.541
e
-
5
2.225
e
-
5
1.900
e
-
5
4.377
e
-
5
SEA-error
1.275
e
-
8
1.442
e
-
8
1.398
e
-
8
1.867
e
-
8Table 2
The condition number.
(
n
1
,
n
2
,
n
3
,
n
4
)
(
2
3
,
2
3
,
2
3
,
2
3
)
(
2
4
,
2
4
,
2
4
,
2
4
)
(
2
5
,
2
5
,
2
5
,
2
5)
|
λ
min
⁡
|
2.214
e
-
3
1.059
e
-
3
5.236
e
-
4
|
λ
max
⁡
|
0.6408
0.6405
0.6404
Cond
2.894
e
+
002
6.047
e
+
002
1.223
e
+
003


## 6. Concluding Remarks

To close this paper, let us make a few concluding remarks.(1)
Evaluation on entries of discrete matrices is very simple and straightforward without any singular integrals by MQM.(2)
The numerical experiments show that MQM retains the optimal convergence orderO(hmax⁡3) and possesses the optimal condition number O(hmin⁡-1) which shows MQM owns the excellent stability. The approximate solutions accuracy order is at least O(hmax⁡5) after splitting extrapolation once, which is a great improvement in accuracy.



---

*Source: 812505-2014-07-10.xml*