The interior field method for Laplace׳s equation in circular domains with circular holes

The interior field method for Laplace׳s equation in circular domains with circular holes

Engineering Analysis with Boundary Elements 67 (2016) 173–185 Contents lists available at ScienceDirect Engineering Analysis with Boundary Elements ...

541KB Sizes 1 Downloads 22 Views

Engineering Analysis with Boundary Elements 67 (2016) 173–185

Contents lists available at ScienceDirect

Engineering Analysis with Boundary Elements journal homepage: www.elsevier.com/locate/enganabound

The interior field method for Laplace's equation in circular domains with circular holes$ Zi-Cai Li a,b, John Y. Chiang b,c, Hung-Tsai Huang d,n, Ming-Gong Lee e,nn a

Department of Applied Mathematics, National Sun Yat-sen University, Kaohsiung 80424, Taiwan Department of Computer Science and Engineering, National Sun Yat-sen University, Kaohsiung 80424, Taiwan c Department of Healthcare Administration and Medical Informatics, Kaohsiung Medical University, 80708, Taiwan d Department of Financial and Computational Mathematics, I-Shou University, Kaohsiung 84001, Taiwan e Department of Leisure and Recreation Management/Ph.D. Program in Engineering Science, Chung Hua University, Hsin-Chu 30012, Taiwan b

art ic l e i nf o

a b s t r a c t

Article history: Received 3 January 2016 Received in revised form 23 February 2016 Accepted 7 March 2016

The null field method (NFM) was proposed by Chen and his co-researchers, and has been discussed in numerous papers, see Chen et al. (2007 [11]), Chen et al. (2002 [12]), Chen et al. (2001 [13]), and Chen and Shen (2009 [14]). The further developments of the NFM have been made in our recent papers (Huang et al., 2013 [21]; Lee et al., 2013 [23]; Lee et al., 2014 [25]; Li et al., 2012 [29]). In this paper, the interior field method (IFM) is proposed, which offers the best performance of the NFM when the field nodes are located exactly on the domain boundary. The algorithms of the IFM are much simpler than those of the NFM, because only one formula of the interior solutions is needed, compared with multiple formulas in the NFM. Since the IFM can be classified into the family of the Trefftz method (Li et al., 2008 [31]), a new error analysis of the IFM and the collocation IFM (CIFM) can be explored, to achieve the optimal convergence rates. Moreover, new proof techniques for aliasing errors in this paper are straightforward, heuristic and much easier to follow, because of direct derivations from trigonometric functions, which are distinct from Canuto and Quarteroni (1982 [8]), Canuto et al. (2006 [9]), and Kreiss and Oliger (1979 [22]). Based on this paper, the IFM and the CIFM may be recommended for those problems solvable by the NFM before. & 2016 Published by Elsevier Ltd.

Keywords: Interior field method Null field method Trefftz method Circular domains Fundamental solutions Error analysis Aliasing errors

1. Introduction For circular domains with circular holes, there exist a number of papers of boundary methods. In Mogilevskaya and Crouch [16,33], the series expansion technique with a direct boundary integral method is used to solve problems involving circular boundary. In Barone and Caulk [5,6] and Caulk [10], the Fourier series are used for the circular holes for boundary integral equations, and in Bird and Steele [7] the simple algorithms as the collocation Trefftz method in [31] are used. In Ang and Kang [2], the complex boundary elements are studied. Recently, the null field method (NFM) has been proposed by Chen and his coresearchers for solving the boundary integral equation in circular ☆ Part of results in this paper were presented at the Sixth TM Workshop jointed with the Second MFS Workshop in National Sun Yat-sen University, Kaoshiung, Taiwan, 15–18 March 2011. n Corresponding author. nn Corresponding author. E-mail addresses: [email protected] (Z.-C. Li), [email protected] (J.Y. Chiang), [email protected] (H.-T. Huang), [email protected] (M.-G. Lee).

http://dx.doi.org/10.1016/j.enganabound.2016.03.006 0955-7997/& 2016 Published by Elsevier Ltd.

domains with circular holes, where the field nodes Q are located exterior to the solution domain S (see (2.18)), and the fundamental solutions (FS) can be expanded as the convergent series. The Fourier series are also used to approximate the known or the unknown Dirichlet and Neumann boundary conditions. Numerous papers have been published for different physical problems by the NFM. Since the algorithms and the error analysis of the IFM are our main concern, we only cite [11–14]. Eq. (2.18) may be called as the exterior field equation, or the null field equation. Based on (2.18) yet without using the Fourier series for the boundary condition, numerical integration by computer is inevitable, to formulate the linear algebraic equations, Tx ¼ b. This approach adversely leads to high computational complexity. Hence, any simplicity of formulation of algebraic equations is significant in applications. Some progress has been made by Wriedt [38], called the null-field method with discrete source (NFM-DS). The transition (T) matrix is provided in Doicu and Wriedt [18], to facilitate the formulation of the algebraic equations. An improving stability of the NFM-DS is reported in [19,38]. The NFM-DS is effective and popular in light scattering. In Li et al. [29], the explicit algebraic equations of the null field method are derived, a strict proof is given for the validity of the

174

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185

field nodes locating on the domain boundary, and a stability analysis is made for simple annular domains. However, even for the Dirichlet problem of Laplace's equation, when the logarithmic capacity (transfinite diameter) C Γ ¼ 1, the solutions may not exist, or not be unique if existing [15], to cause a singularity of the discrete algebraic equations, called the degenerate scales. In Lee et al. [23,25], a comparative analysis is explored for degenerate scales and algorithm singularity of the NFM. The overdetermined system and the truncated singular value decomposition (TSVD) are proposed, to remove the singularity of the discrete matrices obtained, and to achieve good stability of the NFM. Moreover, the NFM has also been applied to Laplace's equation in elliptic domains with elliptic holes in [32,39]. In this paper, the new interior field method (IFM) is proposed. The algorithms of the IFM are much simpler than those of the NFM, because only one formula of the interior solutions is needed, compared with multiple formulas in the NFM. Moreover, since the IFM may be classified into the family of Trefftz methods [31], a new error analysis can be made in this paper, to achieve the optimal convergence rates. The error analysis of the collocation IFM (CIFM) is explored for the simple case of concentric boundaries, where the estimation on aliasing errors is a key issue [9]. New proof techniques for aliasing errors in this paper are straightforward, heuristic and easy for readers to follow. The new proofs based on trigonometric functions are distinct from [8,9,22]. Based on the error analysis and the excellent numerical solutions in this paper, the simple IFM and CIFM may be recommended for those problems solvable by the NFM before (see [21,24,29,32,39]). This paper is organized as follows. In the next section, new algorithms of the IFM are introduced, and a linkage between the IFM and the NFM is explored. In Section 3, new error bounds are derived for the solutions from the IFM, to yield the optimal convergence rates. In Section 4, the error analysis of the collocation IFM (CIFM) is explored for the simple cases of concentric boundaries, and the new proof techniques for aliasing errors are proposed. In the last section, numerical experiments are carried out, and a few concluding remarks are made to address the novelties of this paper.

Fig. 1. The annular domain S.

Fig. 2. The overlapping case of two circular boundaries, ∂SR1 \ ∂SR a ∅.

2. The interior field method 2.1. New algorithms

functions f out and f in in (2.2) are assumed to have approximations of Fourier series,

Consider the Dirichlet problem of Laplace's equation in an annular domain S,

Δu ¼

∂ u ∂ u þ ¼0 ∂x2 ∂y2 2

u ¼ f out

M X

fak cos kθ þbk sin kθg

on ∂SR ;

ð2:3Þ

on ∂SR1 ;

ð2:4Þ

k¼1

2

on ∂SR ;

f out  a0 þ

ð2:1Þ

in S; u ¼ f in

on ∂SR1 ;

ð2:2Þ

where the known functions f out and f in are continuous, and ∂SR and ∂SR1 are the exterior and the interior circular boundaries of S, respectively. The ∂SR with radius R is located at origin ð0; 0Þ, and ∂SR1 with radius R1 is located at origin ðxc ; yc Þ, see Fig. 1. For the Dirichlet problem (2.1) and (2.2), the interior boundary ∂SR1 cannot touch the exterior boundary ∂SR , to guarantee the existence of the solutions. For the overlapping case of the interior boundary, ∂SR1 \ ∂SR a ∅ in Fig. 2, the annular domain S leads to a simply-connected domain, where the minimal interior angle at point P is infinitesimal. When the interior angles of solution domains are infinitesimal, there may not exist any solution of the Dirichlet problems from some counter-examples, see Babuška and Aziz [4] and Li et al. [31, p. 331]. In this paper, we also assume Ra 1, to avoid degenerate scales [25]. Then the solutions of (2.1) and (2.2) exist uniquely. The

f in  a 0 þ

N X

fa k cos kθ þ b k sin kθ g

k¼1

where ak ; bk ; a k and b k are the known coefficients. In (2.3) and (2.4), ðr; θÞ and ðr; θ Þ are the polar coordinates of SR and SR1 with the origins ð0; 0Þ and ðxc ; yc Þ, respectively. Denote nodes x ¼ Q ¼ ðx; yÞ ¼ ðρ; θÞ, and y ¼ P ¼ ðξ; ηÞ ¼ ðr; ϕÞ, wherepxffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ ρ cos θ, y ¼q ρ ffiffiffiffiffiffiffiffiffiffiffiffiffiffi sin θ, ffi ξ ¼ r cos ϕ, and η ¼ r sin ϕ. Then ρ ¼ x2 þ y2 and

r ¼ ξ þ η2 . The fundamental solutions (FS) of Laplace's equation pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi is given by ln j PQ j ¼ ln ρ2  2ρr cos ðθ  ϕÞ þr 2 . The solutions of the Dirichlet problem can be expressed by the boundary integral equation [20], 2

uðxÞ ¼ uðρ; θÞ ¼ 

1 2π

Z  ln j PQ j ∂S

 ∂uðyÞ ∂ ln j PQ j  f ðyÞ dσ y ; ∂ν ∂r

Q A S;

ð2:5Þ

where P A ∂S ¼ ∂SR [ ∂SR1 , and Q is the interior field node inside of S. In (2.5), function f ðyÞ is known from (2.3) and (2.4), but the normal derivatives uν ¼ ∂uðyÞ ∂ν of the solution are unknown. We also

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185

suppose that uν can be approximated by the Fourier series, uν ¼

uν ¼

M  X  ∂u ∂u  p0 þ ¼ pk cos kθ þ qk sin kθ ∂ν ∂r k¼1

on ∂SR ;

N n o X ∂u ∂u ¼   p0 þ p k cos kθ þ q k sin kθ ∂ν ∂r k¼1

on ∂SR1 ;

ð2:6Þ

ð2:7Þ

where pk ; qk ; p k and q k are the unknown coefficients, and ν and ν are the outer normals of ∂SR and ∂SR1 , respectively. Hence, once the coefficients pk ; …; q k are obtained, the final solutions of (2.1) and (2.2) can also be evaluated from integral (2.5), numerically. Under numerical approximation, the direct treatments of (2.5) are troublesome, and will consume a great amount of CPU. We may also employ the series expansions of ln j PQ j (see [1,28]), ln j PQ j ¼ ln j PðyÞ  Q ðxÞj ¼ ln j Pðr; ϕÞ  Q ðρ; θÞj 8 1 X 1 ρ n > i > cos nðθ  ϕÞ; > > U ðx; yÞ ¼ ln r  < n r n¼1

¼ 1 X > 1 r n > > U e ðx; yÞ ¼ ln ρ  cos nðθ  ϕÞ; > : n ρ n¼1

ρ o r;

where the superscripts “e” and “i” designate the exterior and interior field nodes x with respect to ∂SR (or ∂SR1 ), respectively. Based on the orthogonality of trigonometric functions, substituting (2.3), (2.4) and (2.6)–(2.8) into (2.5), we obtain the explicit interior solutions (see [29]), uM  N ¼ uM  N ðρ; θ; ρ ; θ Þ ¼ a0  Rð ln RÞp0 R1 ðln ρ Þp 0

uM  N ðρ; θ ; R1 ; θ Þ ¼ a 0 þ

fa k cos kθ þ b k sin kθ g

on ∂SR1 ;

ð2:13Þ

k¼1

where uM  N ðρ; θ; ρ ; θ Þ is given by (2.9). The algorithms using the interior solutions (2.9) satisfying (2.12) and (2.13) are called the interior field method (IFM) in this paper. We may also establish the collocation equations of (2.12) and (2.13). Since the 2ðM þ N þ 1Þ unknown coefficients in (2.9) are sought, we choose the following 2ðM þN þ 1Þ nodes,

¼ 2M2πþ 1

where Δθ are obtained as

i ¼ 0; 1; …; 2M;

on ∂SR ;

j ¼ 0; 1; …; 2N; and Δθ ¼

2π 2N þ 1.

ð2:14Þ

on ∂SR1 ;

ð2:15Þ

Then the collocation equations

( ) M  pffiffiffiffiffiffi pffiffiffiffiffiffi  Δθ u^ M−N R; θi ; ρi ; θi ¼ Δθ a0 þ ∑ ak cos kθi þ bk sin kθi ; k¼1

i ¼ 0; 1; …; 2M; ( ) o pffiffiffiffiffiffi N n pffiffiffiffiffiffi Δθ u^ M−N ρj ; θj ; R1 ; θj ¼ Δθ a0 þ ∑ ak cos kθj þ bk sin kθj ; j ¼ 0; 1; …; 2N:

ð2:16Þ

Eqs. (2.16) are denoted by the linear algebraic equations,

M  k 1X ρ ðak cos kθ þ bk sin kθÞ þ 2k¼1 R

Tx ¼ b;

þ



N R1 X 1 R1 k ðp k cos kθ þ q k sin kθ Þ 2 k¼1k ρ

þ

N 1X R1 k ða k cos kθ þ b k sin kθ Þ; 2k¼1 ρ

ð2:17Þ nn

ðρ; θÞ A S;

ð2:9Þ

where the transformations of coordinates between ðρ; θÞ and ðρ ; θ Þ are given by qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ρ cos θ−xc ; ð2:10Þ ρ ¼ ðρ cos θ−xc Þ2 þ ðρ sin θ−yc Þ2 ; cos θ ¼ ρ cos θ ¼

k¼1 N X

k¼1

M RX 1ρk þ ðpk cos kθ þ qk sin kθÞ 2k¼1k R

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðρ cos θ þ xc Þ2 þ ðρ sin θ þ yc Þ2 ;

Based on Lemma 2.1, the interior solutions (2.9) are valid until the domain boundary ∂S. Hence the unknown coefficients, pk ; qk ; p k ; q k , can be determined directly from the Dirichlet boundary conditions (2.3) and (2.4), M X fak cos kθ þ bk sin kθg on ∂SR ; ð2:12Þ uM  N ðR; θ; ρ ; θ Þ ¼ a0 þ

ðR1 ; θ j Þ ¼ ðR1 ; jΔθ Þ; ð2:8Þ

ρ¼

Lemma 2.1. Let ðu A H p ð∂SR ÞÞ 4 ðuν A H p  1 ð∂SR ÞÞðp Z 2Þ and ðu A H σ ð∂SR1 ÞÞ 4 ðuν A H σ  1 ð∂SR1 ÞÞðσ Z 2Þ be given. The interior solutions (2.9) hold for the domain boundary ∂S ¼ ∂SR [ ∂SR1 .

ðR; θi Þ ¼ ðR; iΔθÞ;

ρ 4 r;

175

ρ cos θ þ xc : ρ ð2:11Þ

We can prove the following lemma from [24,29].

with n ¼ 2ðM þ N þ 1Þ, the where the discrete matrix T A R unknown vector x A R n consists of the unknown coefficients pk ; qk ; p k ; p k , and the known vector b ∈ Rn . Once the unknown coefficients pk ; qk ; p k ; p k have been found from (2.17), the interior solutions (2.9) become the semi-analytic solutions in S. The algorithms using (2.9) and (2.16) are called the collocation IFM (CIFM) in this paper (Ref. [31]). Note that there are no numerical approximations, and the coefficient matrix T in (2.17) is formulated from (2.16) straightforwardly. The algorithms of the CIFM are very simple and effective. The semi-analytic solutions (2.9) have several advantages over the numerical solutions by other methods, such as the BEM, the NFM-DS, the FEM, the FDM, and the boundary integral methods (BIM): 1. The solutions at any node in S are given explicitly in (2.9). 2. The k-order derivatives at any node in S along any direction can be evaluated from (2.9) easily.

Table 1 Errors and condition numbers for the Model problem by the CIFM. (M,N) J u  u^ M  N J 1;S J u  u^ M  N J 0;∂S J uν  ðu^ M  N Þν J 1;∂S J uν  ðu^ M  N Þν J 0;∂S Cond Cond_eff

(10, 5) 1.39(  4) 1.20(  3) 7.03(  3) 6.61(  3) 1.29(2) 1.06

(20, 10) 1.37(  6) 8.95(  7) 1.20(  5) 9.00(  6) 3.28(3) 2.67

(30, 15) 1.34(  9) 7.28(  10) 1.69(  8) 1.08(  8) 5.86(2) 4.75

(40, 20) 1.30(  12) 6.21(  13) 2.15(  11) 1.22(  11) 8.89(2) 7.19

(50, 25) 1.27(  15) 5.45(  16) 2.59(  14) 1.33( 14) 1.23(3) 9.94

(60, 30) 1.23(  18) 4.88(  19) 2.99(  17) 1.43(  17) 1.61(3) 13.0

176

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185

3. Some engineering analysis, such as flow analysis, stress analysis, profiles of solutions and derivatives in S, etc., can be conducted from (2.9) easily. 4. The semi-analytic solutions (2.9) are most accurate. Take ðM; NÞ ¼ ð60; 30Þ in Model problem in Section 5 for example. Totally, the 92 numerical coefficients are provided in Table 2. The errors and the effective condition numbers are given in Table 1 as ‖u  ðu^ M  N Þ‖1;S ¼ 1:23ð  18Þ; ‖uν  ðu^ M  N Þν ‖1;∂S ¼ 2:99ð  17Þ; Cond_eff ¼ 13:0; to show the high accuracy and the excellent stability. Some stability analysis can be found from [24,29], and the errors analysis of the IFM and the CIFM is first explored in this paper, see Sections 3 and 4.

Proof. Based on Lemma 2.2, it is sufficient to show that (2.12) and (2.13) are equivalent to (2.19) at ∂SR and (2.20) at ∂SR1 , respectively. When ρ ¼ R, Eq. (2.19) leads to Lext ðR; θ; ρ ; θ Þ ¼  Rðln RÞp0 R1 ðlnρ Þp 0 M RX 1 ðp cos kθ þ qk sin kθÞ þ 2k¼1k k

N R1 X 1 R1 k þ ðp k cos kθ þ q k sin kθ Þ 2 k¼1k ρ M 1X ða cos kθ þ bk sin kθÞ 2k¼1 k

N 1X R1 k þ ða k cos kθ þ b k sin kθ Þ ¼ 0: 2k¼1 ρ



2.2. Linkage to the null field method The NFM is formulated from the exterior field equation (or the null field equation) in [20],  Z  ∂uðyÞ ∂ ln j PQ j c ln j PQ j  f ðyÞ ð2:18Þ dσ y ¼ 0; Q A S ; ∂ν ∂ν ∂S where PðyÞ A ∂S, and Q are the exterior field nodes outside of the solution domain S. The explicit algebraic solutions can also be obtained from the orthogonality of trigonometric functions. For the exterior field nodes x ¼ ðρ; θÞ with ρ 4 r ¼ R, the first explicit algebraic equations of the NFM are derived in [29] as Lext ðρ; θ; ρ ; θ Þ ¼  Rðln ρÞp0  R1 ðlnρ Þp 0

M RX 1 R k þ ðpk cos kθ þ qk sin kθÞ 2k¼1k ρ

N R1 X 1 R1 k ðp k cos kθ þ q k sin kθ Þ þ 2 k¼1k ρ M k 1X R ðak cos kθ þ bk sin kθÞ  2k¼1 ρ

N 1X R1 k ða k cos kθ þ b k sin kθ Þ ¼ 0; þ 2k¼1 ρ

ð2:19Þ

þ

ð2:20Þ

After the coefficients pk ; …; q k have been sought from (2.19) and (2.20), the interior solutions (2.9) are also used, to give the final harmonic solutions in S. The algorithms using (2.19), (2.20) and (2.9) are called the null field method (NFM). We can also prove the following lemma from [24,29]. p1

M RX 1 ðp cos kθ þ qk sin kθÞ 2k¼1k k

N R1 X 1 R1 k þ ðp k cos kθ þ q k sin kθ Þ 2 k¼1k ρ

 Rðln RÞp0  R1 ðln ρ Þp 0 þ

M   1X a cos kθ þbk sin kθ 2k¼1 k

N 1X R1 k þ ða k cos kθ þ b k sin kθ Þ ¼ 0; 2k¼1 ρ



ð2:22Þ

Evidently, Eq. (2.21) is identical to (2.22). Similarly, we can show that Lint ðρ; θ; R1 ; θ Þ ¼ 0 is identical to (2.13). This completes the proof of Proposition 2.1.□

Lint ðρ; θ; ρ ; θ Þ ¼  Rðln RÞp0 R1 ðln R1 Þp 0 M RX 1ρk þ ðpk cos kθ þqk sin kθÞ 2k¼1k R

N R1 X 1 ρ k ðp k cos kθ þ q k sin kθ Þ þ a0  a 0 þ 2 k ¼ 1 k R1 M  k 1X ρ ðak cos kθ þ bk sin kθÞ 2k¼1 R

N 1X ρ k ða k cos kθ þ b k sin kθ Þ ¼ 0:  2 k ¼ 1 R1

ð2:21Þ

On the other hand, Eq. (2.12) leads to

ðR; θÞ A ∂SR :

where a common factor 2π is deleted. Next, for the exterior field nodes x ¼ ðρ ; θ Þ with ρ o r ¼ R1 , the second explicit algebraic equations of the NFM are derived in [29] as

p

Proposition 2.1. Let ðu A H p ð∂SR ÞÞ 4 ðuν A H p  1 ð∂SR ÞÞðp Z 2Þ and ðu A H σ ð∂SR1 ÞÞ 4 ðuν A H σ  1 ð∂SR1 ÞÞðσ Z 2Þ be given. The IFM is equivalent to the special NFM when the field nodes Q A ð∂SR [ ∂SR1 Þ.

The first advantage of the IFM over the NFM is simplicity of algorithms, because three different formulas (2.19), (2.20) and (2.9) are employed in the NFM. Moreover, the IFM yields the best numerical performance of the NFM in accuracy and stability, see [21,29]. For the Neumann problem, the IFM is also developed in Lee et al. [24], to also offer the best performance of the second kind of NFM. Since the interior solutions (2.9) are derived directly from the interior field equation (2.5), the interior field method (IFM) is called in this paper, to differentiate from the null field method (NFM). By following (2.5), other numerical treatments are also developed in Ochmann [34] for acoustic radiation and scattering, and the full field method (FFM) is called, to compare to the NFMDS in [18,19,38]. The “interior field” and the “full field” focus on the position of field nodes and the consequence of field effects, respectively. The “explicit interior solutions” (2.9) are essential to both algorithms and analysis. It is due to the “interior” solutions (2.9) that our numerical algorithms (see (2.16)) are very simple, without any numerical integration. Hence, “the interior field method (IFM) and the collocation interior field method (CIFM)” may better represent the characteristics of our numerical methods.

3. New error analysis

σ

Lemma 2.2. Let ðu A H ð∂SR ÞÞ 4 ðuν A H ð∂SR ÞÞðp Z 2Þ and ðu A H ð∂SR1 ÞÞ 4 ðuν A H σ  1 ð∂SR1 ÞÞðσ Z 2Þ be given. Eqs. (2.19) and (2.20) hold for ðρ; θÞ A ∂SR and ðρ ; θ Þ A ∂SR1 , respectively. The intrinsic relation between IFM and NFM is given in the following proposition.

Recently, our efforts are devoted to establish effective numerical algorithms for Laplace's equation in circular and elliptic domains. Several criteria are often used to evaluate the numerical algorithms: (1) accuracy (error bounds and convergence rates), (2) stability (including the algorithm singularity), (3) computation complexity

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185

(CPU time and computer storage), and (4) feasibility of algorithms (including the explicit algebraic equations) and their programming implementation. Among them, accuracy is the most important. The key idea of error analysis is to classify the IFM into the Trefftz family [31]. By following [31], we may derive the error bounds of the IFM. To facilitate the error analysis is another remarkable advantage of the IFM.

For the Dirichlet problem, the interior (i.e., particular) solutions (2.9) are chosen as the admissible functions, where only coefficients, pk, p k , qk and q k , are unknown to be sought. Denote the set of solutions (2.9) by V M  N with ð2ðM þ NÞ þ 2Þ unknown coefficients, and the domain boundary by Γ ¼ ∂SR [ ∂SR1 . Since uM  N satisfies Laplace's equation in S already, the unknown coefficients, pk ; qk ; p k ; q k , can be sought by satisfying the Dirichlet boundary condition only: on Γ ;

ð3:1Þ

where 8 1 X > > > ðak cos kθ þ bk sin kθÞ; > f out ¼ u ¼ a0 þ > < k¼1 f 1 X > > > f in ¼ u ¼ a 0 þ ða k cos kθ þb k sin kθ Þ; > > : k¼1

on

u2 ¼ u2 ðρ ; θ Þ ¼ R1 ðln ρ Þp 0 þ þ

1 1X R1 k in ða k cos kθ þ b k sin kθ Þ≔P in N u2 þ RN u2 : 2k¼1 ρ ð3:10Þ

P M u1 ¼ a0  Rðln RÞp0 þ þ

M  k 1X ρ ðak cos kθ þbk sin kθÞ; 2k¼1 R

ð3:11Þ

1 R X 1ρk ðpk cos kθ þ qk sin kθÞ 2 k ¼ Mþ1 k R 1 ρk 1 X ðak cos kθ þ bk sin kθÞ; þ 2 k ¼ Mþ1 R

P in N u2 ¼  R1 ðln ρ Þp 0 þ

∂SR ;

M RX 1ρk ðpk cos kθ þ qk sin kθÞ 2k¼1k R

RM u1 ¼

ð3:2Þ on



1 R1 X 1 R1 k ðp k cos kθ þ q k sin kθ Þ 2 k¼1k ρ

In (3.9) and (3.10), the notations of the truncated parts and remainders are

3.1. Classification of the IFM into the Trefftz family [31]

u¼f

177

þ

∂SR1 ;

ð3:12Þ



N R1 X 1 R1 k ðp k cos kθ þ q k sin kθ Þ 2 k¼1k ρ

N 1X R1 k ða k cos kθ þ b k sin kθ Þ; 2k¼1 ρ

ð3:13Þ



1 R1 X 1 R1 k ðp k cos kθ þ q k sin kθ Þ 2 k ¼ Nþ1 k ρ k 1 1 X R1 þ ða k cos kθ þb k sin kθ Þ: 2 k ¼ Nþ1 ρ

and ak ; bk ; a k ; b k are the true Fourier coefficients. We will follow the Trefftz method (TM) in [26,31] to seek uM  N . Denote the energy Z ð3:3Þ IðuÞ ¼ ðu  f Þ2 :

Rin N u2 ¼

Then the solution uM  N of the IFM can be obtained by

The smooth solution on ∂S can be expressed by the Fourier series,

Γ

IðuM  N Þ ¼ min IðvÞ: v A VM  N

ð3:4Þ

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi R 2 Denote the zero norm ‖v‖0;Γ ¼ Γ v dℓ. Then the IFM solution uM  N also satisfies ‖u uM  N ‖0;Γ ¼ min ‖u  v‖0;Γ : v A VM  N

ð3:5Þ

When the integral in (3.3) involves numerical approximation, the modified energy is Z c bIðvÞ ¼ ðv  f Þ2 ; ð3:6Þ Γ

R R where b Γ is the numerical approximation of Γ by some quadrature rules, such as the central or the trapezoidal rule. The numerical solution u^ M  N is obtained by bIðu^ M  N Þ ¼ min bIðvÞ: v A VM  N



uj∂SR ¼ f ðθÞ ¼ a○k þ ∑

k¼1

π

0

Denote the truncated Fourier series of f ðθÞ in (3.15) (or the projection solution called) by

k¼1

ð3:17Þ

The Sobolev norms for (3.15) are defined by (Z ‖f ‖0;½0;2π  ¼ j f j 0;½0;2π  ¼

)12



2

f dθ

¼

( )1 i 2 ∞ h pffiffiffi ○ π 2ða○0 Þ2 þ ∑ ða○k Þ2 þ ðbk Þ2 ; k¼1

0

ð3:18Þ 2π

j f j p;½0;2π  ¼



0

p

d f dθp

)12

2

¼



pffiffiffi π

(



2p

∑ k

k¼1

h



ða○k Þ2 þ ðbk Þ2

i

)12 ;

p Z 1;

ð3:19Þ ð3:8Þ

where

( ‖f ‖p;½0;2π  ¼

1 RX 1ρk u1 ¼ u1 ðρ; θÞ ¼ a0  Rðln RÞp0 þ ðpk cos kθ þ qk sin kθÞ 2k¼1k R 1  k 1X ρ þ ðak cos kθ þ bk sin kθÞ≔P M u1 þ RM u1 ; 2k¼1 R

○ ○ ak cos kθ þ bk sin kθ :

M

P M f ¼ P M f ðθÞ ¼ a○0 þ ∑

For the interior solutions in (2.9), we may split it by u ¼ u1 þ u2 ;

ð3:15Þ



(Z

3.2. Preliminary lemmas

○ ○ ak cos kθ þ bk sin kθ ;

where a○k and bk are the true Fourier coefficients defined by Z 2π Z 1 1 2π f ðθÞ dθ; a○k ¼ f ðθÞ cos kθ dθ; a○0 ¼ 2π 0 π 0 Z 2π 1 ○ bk ¼ f ðθÞ sin kθ dθ; k ¼ 1; 2; … ð3:16Þ

ð3:7Þ

Following [31], the collocation equations (2.16) are just equivalent to (3.6), where the trapezoidal rule is used (also see Section 4.1).

ð3:14Þ

p X

k¼0

)12 jf j2k;½0;2π 

:

ð3:20Þ

We have a lemma from [8,9]. ð3:9Þ

Lemma 3.1. Let the solutions (3.15) be given. Suppose f A H p ð0; 2π Þðp Z 1Þ, then there exist the bounds of the remainders

178

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185

of (3.17) ‖f −P M f ‖q;½0;2π  ¼ ‖RM f ‖q;½0;2π

∞ C

○ ¼ ∑ a○ cos kθ þbk sin kθ

q;½0;2π r M p−q j f j p;½0;2π  ;

k ¼ M þ 1 k

0 rq o p;

ð3:21Þ where C is a constant independent of

Theorem 3.1. Let ðu A H p ð∂SR ÞÞ 4 ðuν A H p  1 ð∂SR ÞÞðp Z 2Þ and ðu A H σ ð∂SR1 ÞÞ 4 ðuν A H σ  1 ð∂SR1 ÞÞðσ Z 2Þ be given. When the exact Fourier coefficients of (2.3) and (2.4) are given, the numerical solution from the IFM has the following bound:  ‖u  uM  N ‖0;Γ r C

 1 1 ð‖u‖p;∂SR þ ‖uν ‖p  1;∂SR Þ þ σ ð‖u‖σ ;∂SR þ ‖uν ‖σ  1;∂SR Þ ; 1 1 N Mp

ð3:31Þ

θ and M.

Lemma 3.2. For the remainders (3.12) and (3.14), there exist the bounds,

where u is given in (3.8), ∂SR1 , respectively.

ν and ν are the exterior normals of ∂SR and

‖RM u1 ‖0;∂SR r C‖RM u1 ‖0;∂SR ;

ð3:22Þ

in Proof. Choose v ¼ u M  N ¼ P M u1 þ P in N u2 , where P M u1 and P N u2 are given in (3.11) and (3.13), respectively. From (3.5) we have

in ‖Rin N u2 ‖0;∂SR r C‖RN u2 ‖0;∂SR ;

ð3:23Þ

‖u  uM  N ‖0;Γ r‖u  u M  N ‖0;Γ ¼ ‖u  ðP M u1 þ P in N u2 Þ‖0;Γ

1

1

in ¼ ‖RM u1 þ Rin N u2 ‖0;Γ r‖RM u1 ‖0;Γ þ ‖RN u2 ‖0;Γ ;

where u1 and u2 are given in (3.9) and (3.10). Proof. Define the semi-norm ( )12 Z Z ðvðPÞ  vðQ ÞÞ2 2 dsðPÞdsðQ Þ ; ð3:24Þ ‖v‖1;S ¼ ‖v‖0;S þ 2 ðP  Q Þ2 S S qffiffiffiffiffiffiffiffiffiffiffiffiffiffi R 2 where ‖v‖0;S ¼ S v ds. Since RM u1 in (3.12) are harmonic functions in SR and SR1 ð  SR Þ, we have from [4] ‖RM u1 ‖0;∂SR r C‖RM u1 ‖1;SR r C‖RM u1 ‖1;SR r C 1 ‖RM u1 ‖0;∂SR : 2

1

2

1

ð3:25Þ

1 1

Ω ¼ lim Ωðρ hug Þ;

‖v‖1;Ω ¼ lim ‖v‖1;Ωðρ hug Þ :

ρ hug -1

2

1 1

ρ hug -1

ð3:26Þ

2

Since Rin N u2 in (3.14) are harmonic functions in Ωðρ hug Þ, there exist the bounds from [4], 



in in in ‖Rin N u2 ‖1;Ωðρ hug Þ r C‖RN u2 ‖0;∂Ωðρ hug Þ r C ‖RN u2 ‖0;∂SR þ‖RN u2 ‖0;∂Sρ 2

1

  Since Rin N u2 



ρ hug

‖Rin N u2 ‖0;∂Sρ

hug

¼ O ρ1 hug 0



hug

:

ð3:27Þ

from (3.14) and N Z 0, we have 1

B 1 C ¼ O@qffiffiffiffiffiffiffiffiffiffiA-0;

ρ hug

ρ hug -1:

ð3:28Þ

(

)

in in 1 ‖Rin N u2 ‖1;Ω1 r C ‖RN u2 ‖0;∂SR þ lim ‖RN u2 ‖0;∂Sρ 1

ρ hug -1

hug

ð3:29Þ

¼ C‖Rin N u2 ‖0;∂SR : 1

Also since Rin N u2 are harmonic functions, we obtain from (3.29), in in in 1 1 ‖Rin N u2 ‖0;∂SR r C‖RN u2 ‖1;Ω r C‖RN u2 ‖1;Ω1 r C 1 ‖RN u2 ‖0;∂SR : 2

‖RM u1 ‖0;Γ r ‖RM u1 ‖0;∂SR þ ‖RM u1 ‖0;∂SR r C‖RM u1 ‖0;∂SR 1 8

1

< X 1

ðpk cos kθ þ qk sin kθÞ rC

: k k ¼ M þ1 0;∂S 9 R

X

1 =

ðak cos kθ þ ak sin kθÞ : þ

;

ð3:33Þ

0;∂SR

We obtain from Lemma 3.1 and the orthogonality of trigonometric functions,

( )12

X

1 1 X 1 1 2

2 ðp cos k θ þq sin k θ Þ ¼ π R ðp þ q Þ

k k k 2

k k k ¼ M þ1 k ¼ M þ1 k 0;∂S R

(

)12

1 X 1 πR ðp2k þ q2k Þ M k ¼ Mþ1

1

1

X

¼ ðpk cos kθ þqk sin kθÞ

M

r

k ¼ Mþ1

0;∂SR

1 r C p j uν j p  1;∂SR : M

ð3:34Þ

Also from Lemma 3.1 similarly,

1

X

ðak cos kθ þ bk sin kθÞ

k ¼ M þ1

rC 0;∂SR

1 j uj p;∂SR : Mp

ð3:35Þ

Combining (3.33)–(3.35) gives

From (3.26)–(3.28), we then obtain 2

where RM u1 and are given in (3.12) and (3.14), respectively. By noting Γ ¼ ∂SR [ ∂SR1 , we have from Lemma 3.2

k ¼ M þ1

This is (3.22). In (3.25) and the following text, C and C1 are two bounded constants independent of M and N, and their values may be different in different contexts. 1 1 Denote Ω and Ω1 for two unbounded domains outside of ∂SR 1 1 1 1 and ∂SR1 , respectively. Then Ω  Ω1 and S ¼ Ω1 ⧹Ω . Also define a huge annular domain Ωðρ hug Þ with the interior and the exterior boundaries ∂SR1 and ∂Sρ hug , respectively, where ∂Sρ hug is a circle with a very large radius ρ ¼ ρ hug ð⪢1Þ. Then we have

ð3:32Þ

Rin N u2

2

1

ð3:30Þ

This is (3.23), and completes the proof of Lemma 3.2.□ 3.3. Error bounds In this subsection, we derive the error bounds for the IFM in (3.4), and those for the CIFM (2.16) in the next section. We will analyze the errors by two cases of boundary conditions: (I) ak ; bk ; a k and b k in (2.3) and (2.4) are the true Fourier expansion coefficients; (II) otherwise. First for Case I, we have the following theorem.

‖RM u1 ‖0;Γ r C

 1  j uj p;∂SR þ j uν j p  1;∂SR : Mp

ð3:36Þ

Similarly, from Lemmas 3.1 and 3.2, we obtain in in in ‖Rin N u2 ‖0;Γ r‖RN u2 ‖0;∂SR þ ‖RN u2 ‖0;∂SR1 rC‖RN u2 ‖0;∂SR1 o 1n r C σ j uj σ ;∂SR1 þ j uν j σ  1;∂SR1 : N

ð3:37Þ

The desired result (3.31) follows (3.32), (3.36) and (3.37), and this completes the proof of Theorem 3.1.□ ~ Next for Case II, consider other coefficients a~ k ; b~ k ; a~ k ; b k used in (2.3) and (2.4). Let us design the other Dirichlet problem as

Δu~ ¼ 0; where

in S;

u~ ¼ f~ ;

on ∂S;

ð3:38Þ

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185

179

8  M  ∞ > > > a~ þ ∑ a~ cos kθ þ b~ k sin kθ þ ∑ ðak cos kθ þ bk sin kθÞ; on ∂SR ; > < 0 k¼1 k k ¼ Mþ1 f~ ≡    N  ∞ > ~ > ~ ~ > ak cos kθ þ bk sin kθ ; on ∂SR1 ; > : a0 þ ∑ ak cos kθ þ bk sin kθ þ ∑ k¼1

k ¼ Nþ1

and ak ; bk ; a k ; b k are the true Fourier coefficients of u A ð∂SR [ ∂SR1 Þ. We have the following lemma. Lemma 3.3. There exists the bound, ~ 0;Γ ¼ Oðj η j Þ; ‖u  u‖

ð3:39Þ

where (

η ¼ π R 2ða0  a~ 0 Þ2 þ 2

þ π R1

N h X

o ~ ða k  a~ k Þ2 þ ðb k  b k Þ2  :

ð3:40Þ

The errors of the projection solution are given in Lemma 3.1. Below, we derive the errors of the collocation solutions from (2.16). First, we cite a lemma from [3,17,28].

~ 0;Γ  ‖u~  u~ M  N ‖0;Γ ¼ j η j  oð1Þ; ð3:41Þ ‖u  u~ M  N ‖0;Γ Z ‖u  u‖  a   where a ¼ oðbÞ denotes b ⪡1. To reduce the errors, j η j must be small from (3.41). This leads to a rule for the IFM that the Dirichlet boundary conditions should be given with the exact Fourier coefficients. Remark 3.1. Note that the error bounds in Theorem 3.1 are optimal. For the NFM, only a preliminary argument of convergence is provided in Huang et al. [21]. Since for the harmonic functions, there exist the bounds, ð3:42Þ

the domain errors have the same bounds as those in Theorem 3.1. Moreover, by means of the techniques in Li [27], we have the optimal convergence rates in S,  1 ð‖u‖p;∂SR þ ‖uν ‖p  1;∂SR Þ ‖u uM  N ‖0;S r C 1 Mp þ 2  1 ð‖u‖ þ ‖u ‖ Þ : ð3:43Þ þ σ ;∂SR1 ν σ  1;∂SR1 1 Nσ þ 2

For the error analysis of the collocation IFM (CIFM), the estimation on aliasing errors is a key issue, which is more advanced and challenging than that given in [31]. 4.1. Preliminary lemmas and theorem For the function f ðθÞ, the Fourier series is given by ðak cos kθ þ bk sin kθÞ ¼ P M f ðθÞ þ RM f ðθÞ;

ð4:1Þ

k¼1

○ where ak ¼ a○k and bk ¼ bk are the true Fourier coefficients given in (3.16), and the project solution and the remainder are denoted as P M f ¼ P M f ðθ Þ ¼ a0 þ

M X

ðak cos kθ þ bk sin kθÞ;

k¼1

Z 2π 0

f ðθ Þ dθ ;

ð4:4Þ

(

ΔL ð cos kθÞ ¼

k ¼ 1; 2; …; 2π ;

if

0;

otherwise:

k ¼ νL;

PL  1

i¼0

f ðihÞ with ð4:5Þ

ν ¼ 1; 2; …;

ð4:6Þ

For (4.2), choose the trapezoidal rule T 2M þ 1 ðf Þ with L ¼ 2M þ 1, Zd 2π 0

2M

f ðθÞ dθ ¼ h ∑ f ðihÞ;

ð4:7Þ

i¼0

where h ¼ 2M2πþ 1. We have the following lemma. Lemma 4.2. For the trapezoidal rule in (4.7), there exists the orthogonality for 1 r k; ℓ r M, Zd 2π 0

dθ ¼ 2π ;

Zd 2π 0

Zd 2π 0

Zd 2π

cos kθ cos ℓθ dθ ¼ πδk;ℓ ; sin kθ sin ℓθ dθ ¼ πδk;ℓ ;

cos kθ sin ℓθ dθ ¼ 0;

ð4:8Þ

ð4:9Þ

where δkk ¼ 1 for k Z 1, and δk;ℓ ¼ 0 for k a ℓ.

4. Error bounds of the CIFM

1 X

0

f ðθÞdθ 

ΔL ð sin kθÞ ¼ 0;

0

f ðθÞ ¼ a0 þ

Zd 2π

where the trapezoidal rule, T L ðf Þ ¼ 0 f ðθÞ dθ ¼ h h ¼ 2Lπ . Then there exist the equalities,

Then the numerical solution u~ M  N from the IFM have the errors,

2 ;S

ΔL ðf Þ ¼

Rd 2π

k¼1

‖v‖0;S r C‖v‖1 r C 1 ‖v‖0;Γ ;

ð4:3Þ

Lemma 4.1. For the periodical continuous function f ¼ f ðθÞ on ½0; 2π , denote the errors

) M h i X 2 2 ~ ~ ðak  a k Þ þ ðbk  b k Þ

2ða 0  a~ 0 Þ2 þ

ðak cos kθ þ bk sin kθÞ:

k ¼ Mþ1

k¼1

(

1 X

RM f ¼ RM f ðθÞ ¼

ð4:2Þ

Proof. There exist the trigonometric formulas,  1 cos ðk ℓÞθ þ cos ðk þ ℓÞθ ; 2  1 cos ðk ℓÞθ  cos ðk þ ℓÞθ ; sin kθ sin ℓθ ¼ 2  1 sin ðk þℓÞθ þ sin ðk  ℓÞθ : sin kθ cos ℓθ ¼ 2 cos kθ cos ℓθ ¼

Since 1 rk; ℓ r M, we have k þ ℓ r 2M o L ¼ 2M þ 1. From Lemma 4.1, the trapezoidal rule offers no errors for the integrands, cos kθ cos ℓθ; sin kθ sin ℓθ and sin kθ cos ℓθ with 1 r k; ℓ r M. Since there exists the orthogonality for trigonometric R 2π functions, with respect to the integral 0 ○dθ, the desired results (4.8) and (4.9) are obtained.□ Lemma 4.2 may be regarded as the orthogonality of the discrete trigonometric functions. Now we expand the function f ðθÞ by

180

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185

the interpolation,

þ

I M f ¼ I M f ðθÞ ¼ a^ 0 þ

M X

ða^ k cos kθ þ b^ k sin kθÞ;

ð4:10Þ

M 1X R1 k ða k cos kθ þb k sin kθÞ; 2k¼1 ρ

ðρ; θÞ A S: ð4:19Þ

k¼1

where the interpolation coefficients a^ i and b^ i are obtained by the following ð2M þ 1Þ collocation equations: M X pffiffiffi pffiffiffi hfa^ 0 þ ða^ k cos kθi þ b^ k sin kθi Þg ¼ hf ðθi Þ;

i ¼ 0; 1; …; 2M;

ð4:11Þ

k¼1

where θi ¼ ih and h ¼ 2M2πþ 1. Denote the set of (4.10) as VM. The collocation (i.e., interpolation) solution from (4.11) may read as: I^0 ðu^ M Þ ¼ min I^ 0 ðvÞ;

ð4:12Þ

v A VM

For the Dirichlet conditions, the conditions (2.3) and (2.4) are given in advance. Based on the analysis in Section 3, the known coefficients, ak, bk, a k and b k , are assumed to be the true Fourier coefficient without errors. Since only the numerical coefficients, pk, qk, p k and q k , may be different, the collocation solutions may be denoted as (ref. (4.19)), b u^ M  M ðρ; θÞ ¼ a0  Rðln RÞp^ 0  R1 ðln ρÞp 0 M R X 1ρk þ ðp^ k cos kθ þ q^ k sin kθÞ 2k¼1k R

where I^0 ðvÞ ¼

Zd 2π 0

ðv  f ðθÞÞ2 dθ;

ð4:13Þ

Rd 2π and the integration 0 is defined in (4.7). From Lemma 4.2, the interpolation coefficients are given explicitly by a^ 0 ¼

Z 2π 1 d f ð θ Þ dθ ; 2π 0

1 b^ k ¼

π

Zd 2π 0

a^ k ¼

Z 2π 1d

π

0

f ðθÞ sin kθ dθ;

f ðθÞ cos kθ dθ ;

k ¼ 1; 2; …;

k ¼ 1; 2; …

ð4:15Þ

b^ k ¼ bk :

ð4:16Þ

Hence, only the coefficients a^ k and ak ðk Z1Þ in (4.10) are different. We have a key theorem, whose proof is deferred in the next subsection. Theorem 4.1. Let f ðθÞ A H p ð0; 2π Þðp Z 2Þ satisfy the periodical boundary conditions, f

ðℓÞ

ð0Þ ¼ f

ðℓÞ

ð2π Þ;

ℓ ¼ 0; 1; …; p  1:

ð4:17Þ

For the interpolant solutions from (4.11), there exists the bound of the aliasing errors, ‖I M f  P M f ‖½0;2π  ¼

8 ( M
¼

0

)2 ða^ k  ak Þ cos kθ



k¼1

( M pffiffiffiffi X

π

k¼1

þ

M 1X R1 k ða k cos kθ þ b k sin kθÞ; 2k¼1 ρ

ðρ; θÞ A S:

b , qb in (4.20) are determined by the folThe coefficients p^ k , q^ k , p k k lowing collocation equations: ( ) M  pffiffiffi pffiffiffi  h u^ M−M ðR; ihÞ ¼ h a0 þ ∑ ak cos kih þbk sin kih ; k¼1

i ¼ 0; 1; …; 2M;

ð4:21Þ

( ) o M n pffiffiffi pffiffiffi h u^ M−M ðR1 ; ihÞ ¼ h a0 þ ∑ ak cos kih þ bk sin kih ; k¼1

i ¼ 0; 1; …; 2M;

ð4:22Þ

where h ¼ Δθ ¼ 2M2πþ 1. Eqs. (4.21) and (4.22) are equivalent to (3.7) with the trapezoidal rule (4.7). On the other hand, the approximation (2.6) is modified to the interpolant expansion, uν ðR; θÞ 

M  X  ∂u^ p^ k cos kθ þ q^ k sin kθ ðR; θÞ ¼ p^ 0 þ ∂ν k¼1

on ∂SR ;

where the interpolation coefficients p^ k and q^ k are given by

;

p^ 0 ¼

rC‖f  P M f ‖½0;2π  :

uM  M ðρ; θÞ ¼ a0  Rðln RÞp0  R1 ðln ρÞp 0 R 1ρk þ ðpk cos kθ þ qk sin kθÞ 2k¼1k R M X



M R1 X 1 R1 k ðp k cos kθ þq k sin kθÞ 2 k¼1k ρ

Z 2π 1 d uν ðR; θÞ dθ; 2π 0

p^ k ¼

ð4:23Þ

Z 2π 1d uν ðR; θÞ cos kθ dθ; k ¼ 1; 2; …; π 0 ð4:24Þ

ð4:18Þ

In the subsection, consider the simple case of M ¼N and the concentric boundaries of ∂SR and ∂SR1 (i.e., ρ ¼ ρ and θ ¼ θ Þ. Then the interior solutions (2.9) are simplified as

M  k 1X ρ ðak cos kθ þ bk sin kθÞ þ 2k¼1 R



M R1 X 1 R1 k b b sin kθÞ ðp k cos kθ þ q k 2 k¼1k ρ

91 =2

4.2. Error bounds

þ

þ

ð4:20Þ

)12 ða^ k  ak Þ2

M  k 1X ρ ðak cos kθ þ bk sin kθÞ 2k¼1 R

ð4:14Þ

Eqs. (4.14) and (4.15) are analogous to those of the true Fourier coefficients. Based on Lemma 4.1, there exist the equalities of coefficients for the interpolant and the projection solutions of f ðθÞ, a^ 0 ¼ a0 ;

þ

q^ k ¼

Z 2π 1d

π

0

uν ðR; θÞ sin kθ dθ;

k ¼ 1; 2; …;

ð4:25Þ

Rd 2π and the trapezoidal rule 0 is given by (4.7). Similarly, the approximation (2.7) is modified as uν ðR1 ; θÞ 

M n o X ∂u^ b sin kθ b cos kθ þ q b þ ðR1 ; θ Þ ¼ p p 0 k k ∂ν k¼1

on ∂SR1 ;

ð4:26Þ

b are b and q where the interpolation coefficients p k k Zd 2π b ¼ 1 p uν ðR1 ; θÞ dθ; 0 2π 0

b ¼1 p k

π

Zd 2π 0

uν ðR1 ; θÞ cos kθ dθ;

k ¼ 1; 2; …;

ð4:27Þ b ¼1 q k

π

Zd 2π 0

uν ðR1 ; θÞ sin kθ dθ;

We have the following lemma.

k ¼ 1; 2; …

ð4:28Þ

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185

Lemma 4.3. For the simple case of M¼ N and ðρ; θÞ ¼ ðρ ; θ Þ, the b , in (4.20) are exactly defined b and q interplant coefficients, p^ k ; q^ k ; p k k in (4.23)–(4.28). Proof. The series expansions (2.8) of the FS may be approximated to 8 M 1ρn i > > > U^ ðx; yÞ ¼ ln r− ∑ cos nðθ−ϕÞ; ρ o r; > < n ¼ 1n r

ln j PQ j  n M 1 r e > > > U^ ðx; yÞ ¼ ln ρ− ∑ cos nðθ−ϕÞ; ρ 4 r; > : n ¼ 1n ρ ð4:29Þ with their derivatives, 8 n

M ρ > ^ i ðx;yÞ > 1 > ∂U ∂r ¼ þ cos nðθ−ϕÞ; ∑ > r nþ1 n¼1 r ∂ ln j PQ j <

 M > ∂U^ e ðx;yÞ ∂r r n−1 > > cos nðθ−ϕÞ; ¼− ∑ > ∂r : ρn n¼1

 1 1 ð‖u‖ þ‖u ‖ Þ þ ð‖u‖ þ ‖u ‖ Þ : ν p;∂S p  1;∂S σ ;∂S ν σ  1;∂S R R R1 R1 Mσ Mp ð4:33Þ

Proof. Let Γ ¼ ∂SR [ ∂SR1 , and denote u M  M ¼ P M u1 þ P in M u2 , where P M u1 and P in M u2 are given in (3.11) and (3.13), respectively. We have ‖u  u^ M  M ‖0;Γ r‖u  u M  M ‖0;Γ þ ‖u M  M  u^ M  M ‖0;Γ ^ r ‖RM u1 ‖0;Γ þ ‖Rin M u2 ‖0;Γ þ ‖u M  M  u M  M ‖0;Γ ;

‖u^ M  M  u M  M ‖0;Γ rC

ð4:30Þ

p^ 0 ¼ p0 ;

b ¼p ; p 0 0

u^ M  M  u M  M ¼

ð4:31Þ

b ðx; yÞ are harmonic functions with the same b ðx; yÞ and U where U singularity of ln j PQ j . Substituting (4.29), (4.30), (2.3), (2.4), (4.23) and (4.26) into (4.31), we obtain 8 ) Z 2π ( M X 1
(

) M X 1ρn cos nðθ  ϕÞ R dϕ  ln R  n R n¼1 ( ) Zd M 2π X  a0 þ ½ak cos kϕ þ bk sin kϕ )

M 1 X ρn þ θ  ϕ Þ R dϕ cos nð R n ¼ 1 Rn þ 1 ( ) Zd M 2π X b b b þ p0þ ½p k cos kϕ þ q k sin kϕ 0

 

k¼1 M X n¼1

Rn1  1

ρn

!

)

There exists the bound of (4.37),

X

M 1 ρk

ðp^ k pk Þ cos kθ ‖u^ M  M  u M  M ‖0;Γ r

k R k¼1 :

0;Γ

ð4:38Þ

0;Γ

0;∂SR

0;∂SR1

)2

Z 2π ( X M 1 R1 k ^ þ ðp k  pk Þ cos kθ R1 dθ k R 0 k¼1 ð4:32Þ

Based on the orthogonality in Lemma 4.2, the solutions (4.20) are b , qb are obtained, where the interpolant coefficients p^ k , q^ k , p k k defined in (4.23)–(4.28). This completes the proof of Lemma 4.3.□



M X 1 R1 k ^ k pk Þ2 þ π R1 ð p ðp^ k  pk Þ2 2 2 R k¼1k k¼1k

2

X

M M X 1 1

2 ^ ^ ðp k pk Þ cos kθ r π ðR þR1 Þ ðp  pk Þ r C : 2 k

k k¼1k k¼1 0;½0;2π  ¼ πR

M X 1

ð4:39Þ

We have the following important theorem. Theorem 4.2. Let ðu A H p ð∂SR ÞÞ 4 ðuν A H p  1 ð∂SR ÞÞðp Z 2Þ and ðu A H σ ð∂SR1 ÞÞ 4 ðuν A H σ  1 ð∂SR1 ÞÞðσ Z 2Þ be given. For the simple case of M ¼N and ðρ; θÞ ¼ ðρ ; θ Þ, the collocation solutions (4.20) from (4.21) and (4.22) of the CIFM have the error bound, ‖u  u^ M  M ‖0;∂SR [ ∂SR1

ð4:37Þ

)2 Z 2π ( X M 1 ðp^  pk Þ cos kθ Rdθ ¼ k k 0 k¼1

)

cos nðθ  ϕÞ R1 dϕ :



M R1 X 1 R1 k b ðp k  p k Þ cos kθ; ðρ; θÞ A S: 2 k¼1k ρ

0;Γ

k¼1

0

ð4:36Þ

2

X

M 1ρk

þ ðp^ k  pk Þ cos kθ

k R k¼1

)

M X 1 R1 n  ln ρ  cos nðθ  ϕÞ R1 dϕ n ρ n¼1 ( ) Zd M 2π X þ a0 þ ½a k cos kϕ þ b k sin kϕ (

b ¼q : q k k

From the orthogonality of trigonometric functions, we have

2

2

X

X

M M 1ρk 1ρk



ðp^ k  pk Þ cos kθ ¼ ðp^ k  pk Þ cos kθ



k k R R k¼1 k¼1



(

q^ k ¼ qk ;



X

M 1 R1 k b

þ ðp k p k Þ cos kθ

k ρ k¼1

k¼1

(

ð4:35Þ

M RX 1ρk ðp^ k  pk Þ cos kθ 2k¼1k R

þ

e

0

p; σ Z2:

Then the aliasing errors are obtained from (4.20) and (4.19),

1

i

 1 1 j uν j p  1;∂SR þ σ j uν j σ  1;∂SR1 ; M Mp

Below we prove (4.35). From Lemma 4.1 and (4.16), there exist the equalities for the coefficients in (4.23)–(4.28),

ρ 4 r:

x A S;

ð4:34Þ

where RM u1 and Rin M u2 are given in (3.12) and (3.14), respectively. The desired result (4.33) follows from (3.36) and (3.37) and the following bound of aliasing errors: 

ρ o r;

When the trapezoidal rule (4.7) with h ¼ 2M2πþ 1 is applied to the interior field equation (2.5), we have 8 9 8 Z b i ðx; yÞ= ^ 1
 rC

181

To apply Theorem 4.2 for the last term in (4.39), we define an auxiliary function wðθÞ ¼ p0 þ

1 X 1 ðp cos kθ þ qk sin kθÞ; k k k¼1

ð4:40Þ

182

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185



X

M 1 c 1

r C ðp k  p k Þ cos kθ rC σ j uν j σ  1;∂SR1 : 2

M k k¼1 0;½0;2π 

ð4:48Þ

Hence, the desired result (4.35) is obtained by (4.38), (4.47) and (4.48). This completes the proof of Theorem 4.2. □ Remark 4.1. Theorem 4.2 is valid for the simple case of ðρ; θÞ ¼ ðρ ; θ Þ. For the eccentric boundaries of ∂SR and ∂SR1 , i.e., ðρ; θÞ a ðρ ; θ Þ, the same error bounds as (3.31) in Theorem 3.1 can be derived by following [31], based on an equivalence between 1

2 Rd 2π two norms: ‖v‖½0;2π  and ‖v‖ ½0;2π  , where ‖v‖ ½0;2π  ¼ f 0 v2 dθg .

4.3. New proof of Theorem 4.1 The new proof techniques for aliasing errors in this subsection are distinct from [8,9,22]. First, we cite a lemma from [3,17,28]. Fig. 3. The location of field nodes for Model problem at ðM; NÞ ¼ ð10; 5Þ.

with the derivatives w 0 ðθ Þ ¼

1 X

ðqk cos kθ pk sin kθÞ:

ð4:41Þ

Lemma 4.4. For the periodical continuous function f ¼ f ðθÞ on ½0; 2π , there exist the errors Z 2π ∞   T N ðf ðxÞ cos ðkxÞÞ− f ðxÞ cos ðkxÞdx ¼ π ∑ a○iNk þ a○iN þ k ; i¼1

0

ð4:49Þ

k¼1

where the trapezoidal rule TN is defined in Lemma 4.1, and a○k are defined in (3.16).

We have from (3.19) j w0 ðθÞj 2p;½0;2π  ¼ π

1 X k¼1

  1∂u2 2p k ðp2k þ q2k Þ ¼   R ∂ρ

; p Z1;

ð4:42Þ

p;∂SR

where the coefficient qk and pk are also given in (2.6). Since p1 ∂u ð∂SR Þ, we conclude that w0 ðθÞ A H p  1 ð0; 2π Þ and ∂ρ A H wðθÞ A H p ð0; 2π Þðp Z 2Þ. To compare with (4.40), the function wðθÞ is rewritten as the Fourier series, wðθÞ ¼ α0 þ

1 X

Lemma 4.5. Let f ðxÞ A H p ð0; 2π Þðp Z1Þ boundary conditions, f

ðℓÞ

ð0Þ ¼ f

ðℓÞ

ð4:43Þ

the

periodical

ℓ ¼ 0; 1; …; p  1:

ð4:50Þ

For the coefficients (3.16), there exist the equalities: p

ðαk cos kθ þ βk sin kθÞ;

ð2π Þ;

satisfy

a○k ¼ ð  1Þ2

ða○k ÞðpÞ k

p

;

ð4:51Þ

for even p;

k¼1

a○k ¼ ð−1Þ

where the Fourier coefficients are given as 1 α k ¼ pk ; k

1 β k ¼ qk : k

ð4:44Þ

The interpolation solution of wM ðθÞ of order M is also denoted by ^ M ðθÞ ¼ α^ 0 þ w

M X

ðα^ k cos kθ þ β^ k sin kθÞ;

ð4:45Þ

k¼1

with α^ k ¼ 1kp^ k and β^ k ¼ 1kq^ k . By applying Theorem 4.1 to function wðθÞ, we have



X

X

M M 1



^ ^ ðp k  pk Þ cos kθ ¼ ðα k  αk Þ cos kθ



k k¼1 k¼1 0;½0;2π 

0;½0;2π 

^ N ðθÞ wN ðθÞ‖0;½0;2π  r C‖w  wN ðθÞ‖0;½0;2π  ¼ ‖w rC

1 1 1 j wðθÞj p;½0;2π  ¼ C p j w0 ðθÞj p  1;½0;2π  rC p j uν j p  1;∂SR ; Mp M M ð4:46Þ

where we have also used Lemma 3.1. Combining (4.39) and (4.46) gives

X

M 1ρk 1

ðp^ k  pk Þ cos kθ r C p j uν j p  1;∂SR : ð4:47Þ

k R M k¼1 0;Γ

Similarly, we have



M

X

1 R1 k c

ððp k  p k Þ cos kθÞ

k ρ k¼1

pþ1 2



ðbk ÞðpÞ p ; k ○

where ða○k ÞðpÞ and ðbk ÞðpÞ are the Fourier coefficients of the pth order ðpÞ derivatives f ðθÞ. Proof. We have from the integration by parts, and from the periodical conditions (4.50), Z 2π π a○k ¼ f ðθÞ cos kθ dθ 0 ( ) Z 2π 1 0 ½f ðθÞ sin kθj 20π  ¼ f ðθÞ sin kθ dθ k 0 Z 2π 1 0 f ðθÞ sin kθ dθ ¼ k 0 Z 1 2π ″ ¼ 2 f ðθÞ cos kθ dθ k 0 Z 2π 1 ð3Þ ¼ 3 f ðθÞ sin kθ dθ k 0 Z 1 2π ð4Þ ¼ 4 f ðθÞ cos kθ dθ: ð4:53Þ k 0 The desired results (4.51) and (4.52) follow from (4.53). □ It is ready to prove Theorem 4.1. From Lemma 4.4, we have IM f  P M f ¼

M X



ða^ k a○k Þ cos kθ ¼ π

k¼1 0;Γ

ð4:52Þ

for odd p;

M X 1 X

fa○iN  k þ a○iN þ k g cos kθ;

k¼1i¼1

ð4:54Þ

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185

183

Table 2 The coefficients pk and p k for the Model problem by the CIFM with ðM; NÞ ¼ ð60; 30Þ. k

pk

k

pk

0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30

0.5770780163555853629332801  0.5770780163555853628963886 0.2885390081777926814178821  0.1442695040888963406701890 0.07213475204444817028970534  0.03606737602222408509241770 0.01803368801111204248643038  0.009016844005556021175870040 0.004508422002778010512855215  0.002254211001389005173487203 0.001127105500694502495848805  0.0005635527503472511490058902 0.0002817763751736254675098481  0.0001408881875868126186501646 0.00007044409379340618608179711  0.00003522204689670296163973020 0.00001761102344835134124693208  8.805511724175522868857169(  6) 4.402755862087605491095156(  6)  2.201377931043638608494491(  6) 1.100688965521646969946144(  6)  5.503444827606429509390254(  7) 2.751722413801327399571965(  7)  1.375861206898694317652299(  7) 6.879306034472957411349586(  8)  3.439653017215144113575639(  8) 1.719826508585417008008709(  8)  8.599132542697329697202923(  9) 4.299566271110704450835844(  9)  2.149783135309186634384652(  9) 1.074891567400222437318127(  9)

31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

 5.374457834375349840578695(  10) 2.687228914479858574113285(  10)  1.343614454450058630972183(  10) 6.718072243531041399962209(  11)  3.359036092225722373866789(  11) 1.679518015752515456959824(  11)  8.397589766953642333672007(  12) 4.198794563462410114756047(  12)  2.099396953511326922731225(  12) 1.049698140330340502294145(  12)  5.248487255344506854437989(  13) 2.624240099312103013174767( 13)  1.312116439245049716182107( 13) 6.560545271650360329182492(  14)  3.280234890875886916592870(  14) 1.640078880301636176703591(  14)  8.200000552154518162222303(  15) 4.099598236761487811256666(  15)  2.049388905712056680872966(  15) 1.024276101186992265077194(  15)  5.117116309616405876035722(  16) 2.554214747618674265096856(  16)  1.272687792136474696304491(16) 6.318544164314837382881272(  17)  3.113808037052595679435140(  17) 1.511138850643660065392109(  17)  7.100572662746273063186526(  18) 3.109146617169360886048791(  18)  1.151081619474320145803881(  18) 2.586055412206235907184694(  19)

k

pk

k

pk

0 1 2 3 4 5 6 7 8 9 10 11 12 3 14 15

 1.442695040888963407322006 0.7213475204444817035609846  0.1803368801111204259012696 0.04508422002778010632164796  0.01127105500694502668242173 0.002817763751736256458774890  0.0007044409379340643071336397 0.0001761102344835158025544843  0.00004402755862087922934050628 0.00001100688965521946561349442  2.751722413805227892561753(  6) 6.879306034508935767280875(  7)  1.719826508631657016640818(  7) 4.299566271530320084956572(  8)  1.074891567934832257361661(  8) 2.687228919271517175929575(  9)

16 17 18 19 20 21 22 23 24 25 26 27 28 29 30

 6.718072304210806486453059(  10) 1.679518069600705358405000(  10)  4.198795242524197672962448(  11) 1.049698737901264828976625(  11)  2.624247614217697315702113(  12) 6.560610911576025631657819(  13)  1.6401612981483780910979329(  13) 4.100313093303886134680107( 14)  1.025173240859096409534831(  14) 2.561936644771091403478418(  15)  6.415354275811859716764790(  16) 1.592845270009908268545692(  16)  4.099881720056604575302988(  17) 9.057372598245342328178038(  18)  3.743032071033984720163576(  18)

where N ¼ 2M þ1. Then there exists the bound, ‖I M f  P M f ‖20;½0;2π  r π3

M X

(

k¼1

1 X

¼π

3

M X k¼1

(

1 X

Then we have

)2

ða○n ÞðpÞ -0; as n-1:

ða○iN  k þa○iN þ k Þ

First, assume

i¼1

j a○iN  k j

)2 ðj a○iN  k j þ j a○iN þ k j Þ

:

ð4:55Þ

i¼1

We only show the proof for the even integer p, because the proof for the odd p is similar. From Lemma 4.5 we have 1 X i¼1

j a○iN  k j ¼

1 X

j ða○iN  k ÞðpÞ j p : ðiN kÞ

i¼1

From Parseval's identity and the assumption 2π Þðp Z 2Þ, the infinite series converges, Z 1 X 2 1 2π ðpÞ 2 fða○n ÞðpÞ g ¼ ff ðθÞg dθ r C; n¼1

π

0

r C j a○N  k j ;

f ðθÞ A H p ð0;

there must exist the bounds from (4.58)

i ¼ 1; 2; …;

ð4:59Þ

where C is a constant independent of i, k and N. Hence Eq. (4.56) leads to 1 X

j a○iN  k j r C

i¼1

ð4:56Þ

ð4:58Þ

a○N  k a 0,

ðpÞ 1 X 1 N  kÞ j ¼ C j ða○N  k ÞðpÞ j : ðiN  kÞp ðiN  kÞp i¼1

1 X j ða○ i¼1

ð4:60Þ

For the infinite series in (4.60), we have 1 X

1 X 1 1 1 ; p¼ pþ ðiN  kÞ ðN  kÞ ðiN  kÞp i¼1 i¼2

ð4:61Þ

where ð4:57Þ

1 X

1 r ðiN  kÞp i¼2

Z 1

1

1 1 1 1 1 dx ¼ r : Nðp  1Þ ðN  kÞp  1 p  1 ðN  kÞp ðNx  kÞp

ð4:62Þ

184

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185

Combining (4.61) and (4.62) gives

1 X 1 1 1 p 1 C r 1 þ ¼ r ; p  1 ðN  kÞp p  1 ðN  kÞp ðN  kÞp ðiN  kÞp i¼1

LIFM int ðρ; θ ; ρ ; θ Þ ¼  Rðln RÞp0  R1 ðlnρ Þp 0 þ ð4:63Þ þ

by noting p Z 2. By applying Lemma 4.5 again, we have from (4.60) and (4.63) 1 X

j a○iN  k j rC j ða○N  k ÞðpÞ j

i¼1

1 ¼ C j a○N  k j : ðN  kÞp

ð4:64Þ

Similarly, we have 1 X

j a○iN þ k j rC j a○N þ k j :

ð4:65Þ

i¼1

From (4.64) and (4.65), Eq. (4.55) leads to 

M

‖I M f −P M f ‖20;½0;2π r C ∑

k¼1

j a○N−k j þ j a○N þ k j n

M

r 2C ∑

k¼1

( ¼ 2C

j a○ð2M þ 1Þ−k j 2 þ j a○ð2M þ 1Þ þ k j 2

2M



k ¼ Mþ1

r 2C

2





k ¼ M þ1

ða○k Þ2 þ

3M þ 1

ða○k Þ2 r 2C

¼ 2C‖RM f ‖20;½0;2π 



k ¼ 2M þ 2

n

k ¼ M þ1



ða○k Þ2 þ ðbk Þ2

o

¼ 2C‖f −P M f ‖20;½0;2π :

ð4:66Þ

This is the desired result (4.18). To complete the proof, consider the case that there just happens some zero coefficient with n a○N  kn ¼ 0, 1 r k r M. We may choose the next subsequent nonzero coefficient, say, a○2N  kn a 0. We then replace a○N  k in (4.64) by a○2N  kn , to also reach the same bound (4.66). This completes the proof of Theorem 4.1. □

5. Numerical experiments and concluding remarks

wi LIFM ext ðR; iΔθ ; ρ i ; θ i Þ ¼ 0; i ¼ 0; 1; …; M;

ð5:6Þ

wi LIFM int ðρi ; θ i ; R1 ; iΔθ Þ ¼ 0; i ¼ 0; 1; …; N;

ð5:7Þ

¼ 2M2πþ 1,

2π 2N þ 1,

Δθ ¼ and the collocation nodes are shown where Δθ in Fig. 3. Based pffiffiffi on the stability analysis in [29], the weights, w0 ¼ 1 and wi ¼ 2 for iZ 1, are used in (5.6) and (5.7). Eqs. (5.6) and (5.7) are denoted by the linear algebraic equations, ð5:8Þ

where T A R with n ¼ M þ N þ 2. By using the Gaussian elimination, the coefficients pk and p k are obtained from (5.8). The traditional condition number is defined as Cond ¼ σσmax , where σ max min and σ min are the maximal and the minimal singular values of matrix T in (5.8), respectively. The effective condition number is ‖b‖ defined in [30] as Cond_eff ¼ σ min ‖x‖, where ‖x‖ is the 2-norm of vector x. Based on (5.2), the boundary and the domain errors can be computed easily, see [29]. A better match between M and N is found as ðM; NÞ ¼ ð2; 1Þ, see Wu [37]. Table 1 lists the errors and condition numbers, and Table 2 gives all coefficients at ðM; NÞ ¼ ð60; 30Þ. More computation results are provided in [37]. From Table 1, we can see the following asymptotes: ‖u  ðu^ M  N Þ‖1;S ¼ Oð0:50M Þ;

Cond ¼ OðMÞ;

5.1. Numerical experiments Consider Model problem in [23,29,35] as u¼1

on ∂SR ;

u ¼ 0 on ∂SR1 ;

ð5:1Þ

where R¼ 2.5 and R1 ¼ 1, and the origins of SR and SR1 are located at the origins ð0; 0Þ and ð  1; 0Þ, respectively. The true solution of (5.1) can be found in [14,35] as ( ) 1 16ρ 2 þ 1 þ 8ρ cos θ ln ; ð5:2Þ uðρ ; θ Þ ¼ 2 ln 2 ρ 2 þ 16 þ 8ρ cos θ where ðρ ; θ Þ are the polar coordinates of SR1 with origin ð  1; 0Þ. By means of symmetry, the interior solutions (2.9) is simplified as uM  N ðρ; θÞ ¼ 1  Rðln RÞp0 R1 ðln ρ Þp 0 þ þ

M RX pk ρk cos kθ 2k¼1 k R



N R1 X p k R1 k cos kθ ; ðρ; θÞ A S: 2 k¼1 k ρ

ð5:3Þ

We obtain two boundary equations directly from (5.3) and (5.1), LIFM ext ðρ; θ ; ρ ; θ Þ ¼  Rðln RÞp0  R1 ðlnρ Þp 0 þ þ

M RX pk cos kθ 2k¼1 k



N R1 X p k R1 k cos kθ ¼ 0; 2 k¼1 k ρ

ð5:5Þ

The coefficients, p0 ; p 0 ; pk ; p k , are unknowns, and the total number of unknowns is ðM þ N þ 2Þ. From (5.4) and (5.5), we choose ðM þ N þ 2Þ collocation equations,

‖uν  ðu^ M  N Þν ‖1;∂S ¼ Oð0:51M Þ;

Δu ¼ 0 in S;

ðρ ; θ Þ A ∂SR1 :

nn

ða○k Þ2



N R1 X pk cos kθ þ 1 ¼ 0; 2 k¼1 k

Tx ¼ b;

o

)



M RX pk ρk cos kθ 2k¼1 k R

ðρ; θÞ A ∂SR ;

ð5:4Þ

‖u  ðu^ M  N Þ‖0;S ¼ Oð0:50M Þ; ‖uν  ðu^ M  N Þν ‖0;∂S ¼ Oð0:51M Þ;

Cond_eff ¼ OðMÞ:

ð5:9Þ ð5:10Þ ð5:11Þ

For the infinitely smooth solution of (5.2) (i.e., p; σ -1), the polynomial convergence rates in (3.31) and (4.33) lead to the exponential convergence rates. Hence, Eqs. (5.9) and (5.10) coincide with the error analysis in Sections 3 and 4. Moreover, the condition numbers and the effective condition numbers in (5.11) are also consistent with the analysis of excellent stability in [24,29]. Note that the entries T ij of matrix T in (5.8) are given explicitly from (5.6) and (5.7). The algorithms of the CIFM are simple, and the computer programming is facile. Based on the error analysis in Sections 4 and 5 and the above numerical results, the convergence is fast. Then, small M and N may satisfy the engineering requirements. For the Gaussian elimination, the CPU time is Oðn3 Þ, where n ¼ M þ N þ 2. Hence, the simplicity of algorithms, facile programming and less CPU time are the remarkable advantages of the CIFM for engineering computation. 5.2. Concluding remarks To close this paper, let us make a few concluding remarks, to display the novelties of this paper. 1. The goals of this paper are twofold: the algorithms of the interior field method (IFM) and the new error analysis. The IFM is much simpler than the NFM, because only the interior field solution (2.9) is needed, in contrast to multiple formulas, (2.19), (2.20) and (2.9), used in the NFM. Compared with our previous

Z.-C. Li et al. / Engineering Analysis with Boundary Elements 67 (2016) 173–185

2.

3.

4.

5.

papers, the key contribution of this paper is the strict proof of optimal convergence for the IFM and the CIFM, which is the most important criterion to evaluate numerical methods. The IFM is equivalent to the NFM, when the field nodes Q are just located on the domain boundary ∂SR [ ∂SR1 . The IFM yields the best stability and optimal convergence rates among all Q chosen in the NFM, based on the analysis and computation in [21,29]. Note that the error analysis in this paper is also valid for the specific NFM when the field nodes Q A ð∂SR [ ∂SR1 Þ. By following the Trefftz method in [31], the error bounds of the IFM solutions are provided in Theorem 3.1, to achieve the optimal convergence rates. The error analysis for the collocation (i.e., the interpolation) solutions of the CIFM is made in Section 4. For the simple case of ðρ; θÞ ¼ ðρ ; θ Þ, the error bounds of the CIFM are provided in Theorem 4.2. In Section 4.3, the proofs of Theorem 4.1 distinct from [8,9,22,36] are new and easier to follow. The error analysis in this paper provides a theoretical foundation for the IFM used in [21,24,32,39]. The excellent stability by the NFM can also be achieved by following [24,29,32]. Therefore, the simple IFM and CIFM may be recommended for those problems solvable by the NFM before.

Acknowledgements Authors are indebted to S.R. Wu for the computation in Section 5.1, and grateful to Professors J.T. Chen and C.S. Chen, and the reviewers for their valuable suggestions on the paper.

References [1] Abramowitz M, Stegun IA. Handbook of mathematical functions with formulas, graphs and mathematical tables. New York: Dover Publications, Inc.; 1964. [2] Ang WT, Kang I. A complex variable boundary element method for elliptic partial differential equations in a multiple-connected region. Int J Comput Math 2000;75:515–25. [3] Atkinson EA. An introduction to numerical analysis. second ed. New York: John Wiley & Sons; 1989. [4] Babuška I, Aziz AK. Survey lectures on the mathematical foundations of the finite element method. In: Aziz AK, editor. The mathematical foundations of the finite element method with applications to partial differential equations. New York: Academic Press; 1972. p. 3–359. [5] Barone MR, Caulk DA. Special boundary integral equations for approximate solution of Laplace's equation in two-dimensional regions with circular holes. Q J Mech Appl Math 1981;34:265–86. [6] Barone MR, Caulk DA. Special boundary integral equations for approximate solution of potential problem in three-dimensional regions with slender cavities of circular cross-section. IMA J Appl Math 1985;35:311–25. [7] Bird MD, Steele CR. A solution procedure for Laplace's equation on multiply connected circular domains. Trans ASME 1992;59:398–404. [8] Canuto C, Quarteroni A. Approximation results for orthogonal polynomial in Sobolev spaces. Math Comput 1982;38:67–86. [9] Canuto C, Hussaini MY, Quarteroni A, Zhang TA. Spectral methods, fundamentals in single domains. Berlin: Springer; 2006. [10] Caulk DA. Analysis of steady heat conduction in regions with circular holes by a special boundary-integral method. IMA J Appl Math 1983;30:231–46. [11] Chen JT, Chen CT, Chen PY, Chen IL. A semi-analytical approach for radiation and scattering problem with circular boundaries. Comput Methods Appl Mech Eng 2007;196:2751–64.

185

[12] Chen JT, Lee CF, Chen JL, Lin JH. An alternative method for degenerate scale problem in boundary element methods for two-dimensional Laplace equation. Eng Anal Bound Elem 2002;26:559–69. [13] Chen JT, Lin JH, Kuo SR, Chui YP. Analytical study and numerical experiments for degenerate scale problems in the boundary element method using degenerate kernel and circulants. Eng Anal Bound Elem 2001;25:819–28. [14] Chen JT, Shen WC. Null-field approach for Laplace problems with circular boundaries using degenerate kernels. Numer Methods PDE 2009;25:63–86. [15] Christiansen S. Integral equations without a unique solution can be made useful for solving some plane harmonic equations. J Inst Math Appl 1975;16 (2):143–59. [16] Crouch SL, Mogilevskaya SG. On the use of Somigliana's formula and Fourier series for elasticity problems with circular boundaries. Int J Numer Methods Eng 2003;58:537–78. [17] Davis PF, Rabinowitz P. Methods of numerical integral integration. second ed.. New York: Academic Press; 1983. [18] Doicu A, Wried T. Calculation of the T matrix in the null-field method with discrete sources. J Opt Soc Am A 1999;16:2539–44. [19] Hellmers J, Schmidt V, Wriedt T. Improving the numerical instability of Tmatrix light scattering calculations for extreme articles shapes using the nullfield method with discrete sources. J Quant Spectrosc Radioact Transf 2011;112:1679–86. [20] Hsia GC, Wendland WL. Boundary integral equations. Berlin: Springer; 2008. [21] Huang HT, Lee MG, Li ZC, Chiang JY. Null field and interior field methods for Laplace's equation in actually punctured disks. Abst Appl Anal 2013;2013:1–15 ID Num. 927873. [22] Kreiss HO, Oliger J. Stability of the Fourier method. SIAM J Numer Anal 1979;16:421–33. [23] Lee MG, Li ZC, Huang HT, Chiang JY. Conservative schemes of the null-field method for Dirichlet problems of Laplace's equation. Eng Anal Bound Elem 2013;37:95–106. [24] Lee MG, Li ZC, Huang HT, Chiang JY. Neumann problems of Laplace's equation in circular domains with circular holes by methods of field equations. Eng Anal Bound Elem 2015;51:156–73. [25] Lee MG, Li ZC, Zhang LP, Huang HT, Chiang JY. Algorithm singularity of the null-field method for Dirichlet problems of Laplace's equation in annular and circular domains. Eng Anal Bound Elem 2014;41:160–72. [26] Li ZC. Combined methods for elliptic equations with singularities, interfaces and infinities. Boston: Kluwer Academic Publishers; 1998. [27] Li ZC. Method of fundamental solutions for annular shaped domains. J Comput Appl Math 2009;228:355–72. [28] Li ZC, Huang J, Huang HT. Stability analysis of method of fundamental solutions for mixed boundary value problems of Laplace's equation. Computing 2010;88:1–29. [29] Li ZC, Huang HT, Liaw CP, Lee MG. The null-field method of Dirichlet problems of Laplace's equation on circular domains with circular holes. Eng Anal Bound Elem 2012;36:477–91. [30] Li ZC, Huang HT, Wei Y, Cheng AH-D. Effective condition number for partial differential equations. second ed.. Beijing, China: Science Press; 2015. [31] Li ZC, Lu TT, Hu HY, Cheng AH-D. Trefftz and collocation methods. Southampton, Boston: WIT Press; 2008. [32] Li ZC, Zhang L, Wei Y, Lee MG, Chiang John Y. Boundary methods for Dirichlet problems of Laplace's equation in elliptic domains with elliptic holes. Eng Anal Bound Elem 2015;61:92–103. [33] Mogilevskaya SG, Crouch SL. A Galerkin boundary integral method for multiple circular elastic inclusions. Int J Numer Methods Eng 2001;52:1069–106. [34] Ochmann M. The full field equations for acoustic radiation and scattering. J Acoust Soc Am 1999;105:2574–84. [35] Palaniappan D. Electrostatics of two intersecting conducting cylinders. Math Comput Model 2002;36:821–30. [36] Pasciak JE. Spectral and pseudo spectral methods for advection equations. Math Comput 1980;35:1081–92. [37] Wu SR. The collocation Trefftz method for Laplace's equation on annular shaped domains, circular and elliptic boundaries [Master thesis]. Department of Applied Mathematics, National Sun Yat-sen University, Kaohsiung, Taiwan; June 2011. [38] Wriedt T. Review of the null-field method with discrete source. J Quant Spectrosc Radioact Transf 2007;106:534–45. [39] Zhang L, Li ZC, Lee MG. Boundary methods for mixed boundary problems of Laplace's equation in elliptic domains with elliptic holes. Eng Anal Bound Elem 2016;63:92–104.