# Mechanical Quadrature Methods and Extrapolation for Solving Nonlinear Problems in Elasticity

Authors: Pan Cheng; Ling Zhang
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

This paper will study the high accuracy numerical solutions for elastic equations with nonlinear boundary value conditions. The equations will be converted into nonlinear boundary integral equations by the potential theory, in which logarithmic singularity and Cauchy singularity are calculated simultaneously. Mechanical quadrature methods (MQMs) are presented to solve the nonlinear equations where the accuracy of the solutions is of three orders. According to the asymptotical compact convergence theory, the errors with odd powers asymptotic expansion are obtained. Following the asymptotic expansion, the accuracy of the solutions can be improved to five orders with the Richardson extrapolation. Some results are shown regarding these approximations for problems by the numerical example.


---

## Body

## 1. Introduction

This paper will describe the isolated elastic equations on a bounded planar regionΩ in the plane with nonlinear boundary value conditions:(1)σij,j=0,inΩ,p¯=-g¯x,u¯+f¯x,onΓ,i,j=1,2,where Ω⊂R2 is a connected domain with a smooth closed curve Γ, the stress tensors are σij, n=(n1,n2) is the unit outward normal vector on Γ, the tractor vector is assumed given p¯=p¯1,p¯2⁡T with p¯i=σi1n1+σi2n2, f¯(x)=f¯1x,f¯2xT and continuous on Γ, and g¯(x,u¯)=g¯1x,u¯,g¯2x,u¯T, where g¯i(x,u¯) is a nonlinear function corresponding to the displaced u¯. Following vector computational rules, the repeated subscripts imply the summation from 1 to 2.The problems are studied in many applications, e.g., the circular ring problems and the instance problems [1, 2], and the bending of prismatic bars [3]. Some methods have been proposed for solving the elasticity. The nonconforming mixed finite element methods are established by Hu and Shi [4] to solve the elasticity with a linear boundary. Talbot and Crampton [5] use a pseudo-spectral method to approach matrix eigenvalue problems which is transformed from the governing partial differential equations. The least-square methods are introduced by Cai and Starke [6] for obtaining the solution of the elastic problems. Chen and Hong [7] solved the hyper singular integral equations applying the dual boundary element methods in elasticity. Kuo et al. [8] solve true and spurious eigensolutions of the circular cavity problems that used dual methods. Li and Nie [9] researched the stressed axis-symmetric rods problems by a high-order integration factor method. And the mechanical quadrature methods (MQMs) are adopted by Cheng et al. [10] to solve Steklov eigensolutions in elasticity to obtain high accuracy solutions.Equations (1) are converted into the boundary integral equations [10–12] (BIEs) with the Cauchy and logarithmic singularities by the variational method.(2)δij2u¯jy+∫Γkij∗y,xu¯jxdsx=∫Γhij∗y,xp¯jxdsx,y=y1,y2∈Γ,where δij is the Kronecker delta for i,j=1,2, and the kernels: (3)hij∗=18πμ1-ν-3-4νδijln⁡r+r·ir·j,kij∗=14π1-νr∂r∂n1-2νδij+2r·ir·j+1-2νnir·j-njr·i.The kernels are Kelvin’s fundamental solutions [11]. The Poisson ratio is ν with ν=λ/[2(λ+μ)], r·i means the derivative to xi, and r=y1-x12+y2-x22 is the distance of x and y. The parts ∫Γkij∗(y,x)u¯j(x)dsx of (2) are the Cauchy singularity and the parts ∫Γhij∗(y,x)u¯j(x)dsx are the logarithmic singularity.The nonlinear integral equations are obtained after the boundary conditions are substituted into (2):(4)δij2u¯jy+∫Γkij∗y,xu¯jxdsx=∫Γhij∗y,xg¯jx,u¯dsx+∫Γhij∗y,xf¯jxdsx.This kind of nonlinear integral equation has been discussed in many papers. Rodriguez [12] established sufficient conditions for the existence of solutions to nonlinear, discrete boundary value problems. Wavelet-Galerkin’s methods constructed from Legendre wavelet functions are used by Maleknejad and Mesgarani [13] to approximate the solution. Atkinson and Chandler [14] presented product Simpson’s rule methods and two-grid methods to solve the nonlinear equations. Pao [15] gave a systematic treatment of a class of nonlinear elliptic differential equations and their applications of these problems. Abels et al. [16] discussed the convergence in the systems of nonlinear boundary conditions with Hölder space.Since this is a nonlinear system including the logarithmic singularity and the Cauchy singularity, the difficulty is to obtain the discrete equations appropriately. The displacement vector and stress tensor inΩ can be calculated [17, 18] as follows after the discrete nonlinear equations are solved:(5)u¯iy=∫Γhij∗y,xp¯jxdsx-∫Γkij∗y,xu¯jxdsx,∀y∈Ω,σijy=∫Γhijl∗y,xp¯lxdsx-∫Γkijl∗y,xu¯lxdsx,∀y∈Ω,where (6)kijl∗=1-2νr.jδli+r.iδlj-r.lδij+2r.ir.jr.l4π1-νr,hijl∗=μ2π1-νr22∂r∂n1-2νr.lδij+νr.jδil+r.iδjl-4r.ir.jr.l+2νnir.jr.l+njr.ir.l+1-2ν2nlr.jr.i+njδil+niδjl-1-4νnlδij.Richardson extrapolation algorithms (EAs) are pretty effective parallel algorithms which are based on asymptotic expansion about errors. The solutions, which are solved on coarse grid and fine grid, are used to construct high accuracy solutions. The method is known to be of good stability and optimal computational complexity. Cheng et al. [19] harnessed extrapolation algorithms to obtain high accuracy order for Steklov eigenvalue of Laplace equations. Huang and Lü established extrapolation algorithms to obtain high accuracy solutions for solving Laplace equations on arcs [20] and plane problem in elasticity [21].The reminder of this paper is stated as follows: In Section2 the MQMs are constructed for the nonlinear boundary integral equations and obtain a discrete nonlinear system. In Section 3 the asymptotically compact convergence is proved. In Section 4 an asymptotic expansion of the error about uh is analyzed and the Richardson extrapolation is applied to achieve the accuracy orders O(h5). In Section 5 some results are shown regarding approximations for problems by the numerical example.


## 2. Mechanical Quadrature Methods

We firstly define the notation of the integral operators on boundaryΓ as follows to simplify the equations:(7)Kijuy=∫Γkij∗y,xuxdsxy∈Γ,i,j=1,2,Hijuy=∫Γhij∗y,xuxdsxy∈Γ,i,j=1,2.So (2) can be simplified as the following operator equations:(8)12I0+K11K12K2112I0+K22u¯1u¯2=H11H12H21H22g¯1+f¯1g¯2+f¯2,where I0 is an identity operator.SupposeC2m[0,2π] be the set of 2m times differentiable periodic functions in which the periods of all functions are 2π. A regular parameter mapping [0,2π]→Γ is introduced, so the smooth closed curve x(s)=(x1(s),x2(s)) satisfying x′s2=x1′s2+x2′s2>0 with xi(s)∈C2m[0,2π], i=1,2.The singularity of the kernelshij∗, kij∗ will be analyzed as t→s before the application of discrete methods. Since it is stated in the page 497 of the paper [21] that(9)nir·j-njr·ir=-1i1+Ot-st-s+Ot-s,i≠j,so the Cauchy singularity of kij∗ comes from the part (nir·j-njr·i)/r. Moreover, the logarithmic singularity of hij∗ comes from the component log⁡x-y and the other parts of the operators are smooth. So the operators Hi,j and Ki,j will be divided into several parts.The logarithmic singular operatorHij will be divided into three parts on C2m[0,2π] as follows: (10)A0ut=∫02πa0t,τuτx′τdτ,with a0(t,τ)=c¯0ln⁡2e-1/2sint-τ/2, c¯0=-(3-4ν)/[8πμ(1-ν)], and (11)B0ut=∫02πb0t,τuτx′τdτ,with b0(t,τ)=c¯0(ln⁡xt-xτ-ln⁡2e-1/2sin⁡t-τ/2), and (12)Bijut=∫02πbijt,τuτx′τdτ,with bij(t,τ)=c1r·ir·j, c1=1/[8πμ(1-ν)].The Cauchy singular operatorKij will also be divided into three parts on C2m[0,2π] as follows: (13)C0ut=∫02πc0t,τuτx′τdτ,with c0(t,τ)=c2(n1r·2-n2r·1)/r,c2=-(1-2ν)/[4π(1-ν)], and (14)Miiut=∫02πmiit,τuτx′τdτ,i=1,2,with mii(t,τ)=c3∂r/∂n[(1-2v)+2r·ir·i]/r,c3=-1/[4π(1-ν)], and (15)Mijut=∫02πmijt,τuτx′τdτ,i,j=1,2,i≠j,with mij(t,τ)=c3∂r/∂n(2r·ir·j)/r.We can find thatB0, Bi,j, Mi,j are smooth operators, A0 is a logarithmic singular operator, and C0 is a Cauchy singular operator.Then the equivalent equation of (8) will be shown as(16)12I+C+Mu=A+Bgu+f,where u(t)=(u¯1(x(t)),u¯2(x(t))), f(t)=f1,f2T with fi=Hi,jf¯j(x(t)), g(u(t))=g¯(x(t),u(t)), and (17)I=I000I0,C=0C0-C00,A=A000A0,B=B0+B11B12B21B0+B22,M=M11M12M21M22.Suppose the boundaryΓ will be divided into 2n,(n∈N) equal parts, so the mesh width will be h=π/n and the nodes will be tj=τj=jh, (j=0,1,…,2n-1). The smooth integral operators with the period 2π can easily obtain high accuracy Nyström’s approximations [22]. For example, the approximation operator B0h of B0 can be approximated as follows:(18)B0hut=h∑j=02n-1b0t,τjuτj,and the error will be(19)B0ut-B0hut=Oh2m.Nyström’s approximation Bijh of Bij and Mijh of Mij can be approximated similarly.In order to approximate the logarithmic singular operatorA0, the continuous approximation kernel ap(t,τ) can be defined as follows:(20)apt,τ=a0t,τ,fort-τ≥h,c¯0hln⁡e-1/2h2π,fort-τ<h.By Sidi’s quadrature rules [3], the approximation operator can be constructed as follows:(21)A0hut=h∑j=02n-1apt,τjuτjx′τj,by which the error estimate will be(22)A0ut-A0hut=2∑μ=1m-1ς′-2μ2μ!u2μth2μ+1+Oh2m,where ς′(t) is the derivative of Riemann zeta function.For the Cauchy singular operatorC0, Lü and Huang [21] and Sidi [22] have proposed the operator C0h to approximate it as shown in Lemma 1.Lemma 1.
The Nyström approximate operatorC0h can be defined as(23)C0huti=2c2a1ti,tjh∑j=02n-1cottj-ti2utjx′tjεij,where a1(t,s)=2t-s+Ot-s-1tan⁡((t-s)/2) when s→t, and (24)εij=1,ifi-jis  odd  number,0,ifi-jis  even  number,so the error estimate will be(25)C0uti-C0huti=Oh2m.Thus (16) can be rewritten as follows:(26)12I+Ch+Mhuh-Ah+Bhguh=fh,where Ah, Bh, Ch, and Mh are the approximate matrixes corresponding to the operators A, B, C, and M, respectively.


## 3. Asymptotically Compact Convergence

The Cauchy approximate operator combined with the identity operator will be studied for the reversibility. A property about the operator(1/2)I+Ch will be presented [10].Lemma 2.
( 1 / 2 ) I + C h is invertible and ((1/2)I+Ch)-1 is uniformly bounded.From Lemma2, we can find that the eigenvalues of both C and Ch do not include -1/2. So (16) and (26) can be rewritten as follows: find u∈V(0) satisfying(27)I+L1u+L2gu=f^,where L1=(I/2+C)-1M, L2=(I/2+C)-1(A+B), f^=(I/2+C)-1f, and the space V(m)=C(m)[0,2π]×C(m)[0,2π], m=0,1,2,….Equation (26) can be rewritten as follows:(28)I+L1huh+L2hguh=f^h,with L1h=(I/2+Ch)-1Mh, L2h=(I/2+Ch)-1(Ah+Bh), and f^h=(I/2+Ch)-1fh.Theorem 3.
The approximate operator sequences{L1h} and {L2h} are asymptotically compact sequences and convergent to L1 and L2 in V(0), respectively; i.e.,(29)L1h⟶a.cL1,L2h⟶a.cL2,where →a.c means the asymptotically compact convergence.Proof.
An asymptotically compact sequence will be proved for operator{L2h:V(0)→V(0)}, firstly.
We can find that the kernels ofB0 and Bij(i,j=1,2,) are continuous functions, so the collectively compact convergence [10, 20] is true: (30)B0h⟶c.cB0,Bijh⟶c.cBijinC0,2π,asn⟶∞.Since we have constructed the continuous approximate function of a(t,τ) of ap(t,τ), the approximate operator {A0h} is the asymptotically compact convergent to A0; i.e., A0h→a.cA0 in C[0,2π], as n→∞.
Then we haveAh→a.cA and Bh→a.cB in V(0). It can be concluded that for any bounded sequence {ym∈V(0)} a convergent subsequence must be found in {(Ah+Bh)ym}. Without loss of generality, we assume (Ah+Bh)ym→z, as m→∞. According to the asymptotically compact convergence and the errors quadrature formulae [17, 21], we obtain (31)L2hym-12I+C-1z≤12I+Ch-1Ah+Bhym-z+12I+Ch-1C-Ch12I+C-1z⟶0,asm⟶∞,h⟶0,where · is the norm of £V0,V0. We proved that L2h:V0→V0 is an asymptotically compact sequence.
Furthermore, we will show that{L2h} will be pointwise convergent to L2, as n→∞. We can see that Ah+Bh→a.cA+B∀y∈V(0), so we obtain(32)Ah+Bhy-A+By⟶0,ash⟶0.
From Lemma2, ((1/2)I+Ch)-1 is uniformly bounded, and the errors of the approximate operators of Ah,Bh,Ch will be substituted into the following equations, then (33)L2hy-L2y≤12I+Ch-1·Ah+Bhy-A+By+12I+Ch-1Ch-C12I+C-1A+By⟶0,ash⟶0.So we proved that {L2h} is pointwise convergent to L2, as n→∞.
Since the series{L2h} are asymptotically compact sequences and pointwise convergent to L2, it has been proved that L2h→a.cL2. Similar proof can be done for {L1h}, so the proof of Theorem 3 is completed.In order to obtain the asymptotically compact convergence of the nonlinear equations, some assumptions should be given [14, 23] firstly.Assumptions 4.
(1):g(·,u) is measurable and differentiable for u∈R, and g(x,·) is continuous for x∈Γ.
(2):∂/∂ug(x,u) is Borel measurable and satisfies the inequality 0<l<∂/∂ug(x,u)<L<∞.The asymptotically compact convergence will be drawn according to the assumption about the functiong(u): (34)Ah+Bhg′u⟶a.cA+Bg′u,I+Ch-1Mh+Ah+Bhg′u⟶a.cI+C-1M+A+Bg′u.


## 4. Asymptotic Expansions and Extrapolation

In this section, we derive the asymptotic expansion of errors for the solution and construct the extrapolation algorithm to obtain higher accuracy order solutions.


### 4.1. Asymptotic Expansions

Theorem 5.
Considering the asymptotic property andx(t),f(t)∈V2m[0,2π], there exists a function ω1¯∈V2m-2[0,2π] independent of h, such that(35)u-uht=tj=h3ω¯1t=tj+Oh5Proof.
Let (16) subtract (26) at t=ti, then(36)I2+C+Mu-A+Bgu-I2+Ch+Mhuh+Ah+Bhguh=f-fh;that is,(37)u-uh2+Cu-Chuh+Mu-Mhuh-A+Bgu-Ah+Bhguh=f-fh.
Consider the componentsCu-Chuh in the equations (38)Cu-Chuh=Cu-Chu+Chu-Chuh=C-Chu+Chu-uh,and similar results can be obtained as Mu-Mhuh=(M-Mh)u+Mh(u-uh). Moreover, using the mean value theorem of differentials, we obtain (39)A+Bgu-Ah+Bhguh=A+B-Ah-Bhgu+Ah+Bhgu-guh=A+B-Ah-Bhgu+Ah+Bhg′uh¯u-uhwhere uh¯=uh+t(u-uh)(0≤t≤1).
Note thatM and B are continuous operators, A is a logarithmic operator, and C is a Cauchy operator, so the approximate errors will be (C-Ch)u=O(h2m), (M-Mh)u=O(h2m), f-fh=h3ω1+O(h5), and (40)A+B-Ah-Bhgu=h3ω2+Oh5att=ti.
Equation (37) can be simplified as(41)12I+Ch+Mh-Ah+Bhg′uh¯u-uh=h3ω+Oh5,where ω=ω1+ω2. Since ((1/2)I+Ch)-1 exists and is uniformly bounded, (41) is rewritten as(42)I+Lhu-uh=h3ϕ+Oh5,where Lh=((1/2)I+Ch)-1(Mh-Ah+Bhg′(uh¯)) and ϕ=((1/2)I+Ch)-1ω.
An auxiliary equation will be introduced with the solutionω1¯(43)I+Lω1¯=ϕ,and the corresponding approximate equations will be(44)I+Lhω¯1h=ϕh.Substituting (44) into (42) yields the equations:(45)I+Lhuh-u-h3ω¯1ht=tj=Oh5.Noticing that ω¯1h∈V2m-2[0,2π] and using the results of (42) yield(46)I+Lhω1¯-ω¯1hti=Oh3.Replace ω¯lh by ω¯l and consider Lh is an asymptotic compact operator,(47)uh-u-h3ω1¯t=tj=Oh5;we complete the proof.


### 4.2. Extrapolation Algorithms

An asymptotic expansion about the errors in (35) implies that the extrapolation algorithms [24] can be applied to the solution of (2) to improve the approximate order. The high accuracy order O(h5) can be obtained by computing some coarse grids and fine grids on Γ in parallel. The EAs are described as follows.Step 1.
Take the mesh widthsh and h/2 to obtain the solution of (26) in parallel, where uh(ti),uh/2(ti) are the solutions on Γ corresponding to h and h/2, respectively.Step 2.
Use the values at coarse grids and fine grids to calculate the approximate values atti(48)uh∗ti=178uh/2ti-uhti,and the error is uh∗ti-uti=O(h5).Step 3.
So we can construct a posteriori error estimate:(49)uti-uh/2ti≤uti-178uh/2ti-uhti+17uh/2ti-uhti=uti-uh∗ti+17uh/2ti-uhti≤17uh/2ti-uhti,which is useful for constructing self-adaptive algorithms.


## 4.1. Asymptotic Expansions

Theorem 5.
Considering the asymptotic property andx(t),f(t)∈V2m[0,2π], there exists a function ω1¯∈V2m-2[0,2π] independent of h, such that(35)u-uht=tj=h3ω¯1t=tj+Oh5Proof.
Let (16) subtract (26) at t=ti, then(36)I2+C+Mu-A+Bgu-I2+Ch+Mhuh+Ah+Bhguh=f-fh;that is,(37)u-uh2+Cu-Chuh+Mu-Mhuh-A+Bgu-Ah+Bhguh=f-fh.
Consider the componentsCu-Chuh in the equations (38)Cu-Chuh=Cu-Chu+Chu-Chuh=C-Chu+Chu-uh,and similar results can be obtained as Mu-Mhuh=(M-Mh)u+Mh(u-uh). Moreover, using the mean value theorem of differentials, we obtain (39)A+Bgu-Ah+Bhguh=A+B-Ah-Bhgu+Ah+Bhgu-guh=A+B-Ah-Bhgu+Ah+Bhg′uh¯u-uhwhere uh¯=uh+t(u-uh)(0≤t≤1).
Note thatM and B are continuous operators, A is a logarithmic operator, and C is a Cauchy operator, so the approximate errors will be (C-Ch)u=O(h2m), (M-Mh)u=O(h2m), f-fh=h3ω1+O(h5), and (40)A+B-Ah-Bhgu=h3ω2+Oh5att=ti.
Equation (37) can be simplified as(41)12I+Ch+Mh-Ah+Bhg′uh¯u-uh=h3ω+Oh5,where ω=ω1+ω2. Since ((1/2)I+Ch)-1 exists and is uniformly bounded, (41) is rewritten as(42)I+Lhu-uh=h3ϕ+Oh5,where Lh=((1/2)I+Ch)-1(Mh-Ah+Bhg′(uh¯)) and ϕ=((1/2)I+Ch)-1ω.
An auxiliary equation will be introduced with the solutionω1¯(43)I+Lω1¯=ϕ,and the corresponding approximate equations will be(44)I+Lhω¯1h=ϕh.Substituting (44) into (42) yields the equations:(45)I+Lhuh-u-h3ω¯1ht=tj=Oh5.Noticing that ω¯1h∈V2m-2[0,2π] and using the results of (42) yield(46)I+Lhω1¯-ω¯1hti=Oh3.Replace ω¯lh by ω¯l and consider Lh is an asymptotic compact operator,(47)uh-u-h3ω1¯t=tj=Oh5;we complete the proof.


## 4.2. Extrapolation Algorithms

An asymptotic expansion about the errors in (35) implies that the extrapolation algorithms [24] can be applied to the solution of (2) to improve the approximate order. The high accuracy order O(h5) can be obtained by computing some coarse grids and fine grids on Γ in parallel. The EAs are described as follows.Step 1.
Take the mesh widthsh and h/2 to obtain the solution of (26) in parallel, where uh(ti),uh/2(ti) are the solutions on Γ corresponding to h and h/2, respectively.Step 2.
Use the values at coarse grids and fine grids to calculate the approximate values atti(48)uh∗ti=178uh/2ti-uhti,and the error is uh∗ti-uti=O(h5).Step 3.
So we can construct a posteriori error estimate:(49)uti-uh/2ti≤uti-178uh/2ti-uhti+17uh/2ti-uhti=uti-uh∗ti+17uh/2ti-uhti≤17uh/2ti-uhti,which is useful for constructing self-adaptive algorithms.


## 5. Numerical Example

Example 1.
We first introduce some notations fori=1,2: eih(P)=uihP-uiP is the error of the displacement; rih(P)=eih(P)/eih/2(P) is the error ratio; e¯ih(P)=uih∗P-ui∗P is the error after one-step EAs; and pih(P)=1/7uih/2(P)-uih(P) is a posteriori error estimate.SupposeΩ is an isotropic elliptical body with the axis a=0.3,b=0.5 in the plane domain shown in Figure 1. The parameter formulae for the boundary Γ will be described as x=0.3cos⁡(t), y=0.5sin⁡(t), t∈[0,2π]. And the nonlinear boundary values p=g1,g2T are set as g1=-u12+cos⁡(t), g2=sin⁡(u2)+2sin⁡(t), t∈[0,2π].Figure 1
The shape of the elliptical body.We firstly calculate the numerical solutionsuh=u1h,u2hT on the boundary Γ following (26). Table 1 lists the approximate values of u1h(P) at points P1=(acos⁡π/8,bsin⁡π/8) and P2=(acos⁡π/4,bsin⁡π/4). Table 2 lists the approximate values of u2h(P) at points P1=(acos⁡π/8,bsin⁡π/8) and P2=(acos⁡π/4,bsin⁡π/4).Table 1
The errors and errors ratio ofu1h(P) at points P=P1, P2.
n 16 32 64 128 256 512 e 1 h ( P 1 ) 8.466 E - 4 1.040 E - 4 1.296 E - 5 1.618 E - 6 2.022 E - 7 2.527 E - 8 r 1 h ( P 1 )  8.137 8.028 8.011 8.002 8.000 e ¯ 1 h ( P 1 )  4.12 E - 07 1.29 E - 08 4.03 E - 10 1.24 E - 11 3.93 E - 13 r ¯ 1 h ( P 1 )   31.92 31.99 32.41 31.68  e 1 h ( P 2 ) 6.085 E - 4 7.501 E - 5 9.292 E - 6 1.160 E - 6 1.450 E - 7 1.813 E - 8 r 1 h ( P 2 )  8.112 8.073 8.008 8.000 8.000 e ¯ 1 h ( P 2 )  3.91 E - 07 1.24 E - 8 4.00 E - 10 1.30 E - 11 4.03 E - 13 r ¯ 1 h ( P 2 )   31.47 31.09 30.85 32.16Table 2
The errors and errors ratio ofu2h(P) at points P=P1, P2.
n 16 32 64 128 256 512 e 2 h ( P 1 ) 4.169 E - 4 5.100 E - 5 6.309 E - 6 7.805 E - 7 9.753 E - 8 1.219 E - 8 r 2 h ( P 1 )  8.175 8.083 8.018 8.003 8.000 e ¯ 2 h ( P 1 )  7.67 E - 06 2.44 E - 07 7.80 E - 9 2.45 E - 10 7.75 E - 12 r ¯ 2 h ( P 1 )   31.44 31.29 31.77 31.65  e 2 h ( P 2 ) 9.338 E - 4 1.147 E - 4 1.420 E - 5 1.774 E - 6 2.217 E - 7 2.771 E - 8 r 2 h ( P 2 )  8.138 8.080 8.005 8.001 8.000 e ¯ 2 h ( P 2 )  8.30 E - 06 2.75 E - 7 8.75 E - 9 2.79 E - 10 8.83 E - 12 r ¯ 2 h ( P 2 )   30.21 31.40 30.69 31.56From Tables1-2, we can numerically see (50)log2rihP≈3,log2r¯ihP≈5,which shows that the convergent orders are three for the approximation solutions uih and will be improved to five orders after EAs.The normal derivativeph=p1h,p2hT on Γ can be obtained when the displacement uh=(u1h,u2h)T on Γ is substituted into the boundary condition of (1). So following (4), the displacement uh=u1h,u2hT in Ω can be calculated. Table 3 shows the approximate values of the displacement u1hP0,u2hP0T at an inner point P0=(1/8)(cos⁡π/3,sin⁡π/3) in Ω.Table 3
The errors and errors ratio ofu1h,u2hTat inner points P=P0.
n 8 16 32 64 128 256 e 1 h ( P 0 ) 8.156 E - 4 9.700 E - 5 1.172 E - 5 1.447 E - 6 1.808 E - 7 2.260 E - 8 r 1 h ( P 0 )  8.409 8.276 8.099 8.002 8.000 e ¯ 1 h ( P 0 )  2.13 E - 05 6.81 E - 07 2.21 E - 08 7.04 E - 10 2.22 E - 11 p 1 h ( P 0 )  9.727 E - 5 1.183 E - 5 1.452 E - 6 1.809 E - 7 2.260 E - 8  e 2 h ( P 0 ) 4.372 E - 4 5.145 E - 5 6.227 E - 6 7.684 E - 7 9.599 E - 8 1.200 E - 8 r 2 h ( P 0 )  8.498 8.262 8.104 8.005 8.000 e ¯ 2 h ( P 0 )  4.82 E - 05 1.58 E - 6 5.09 E - 8 1.62 E - 9 5.11 E - 11 p 2 h ( P 0 )  5.154 E - 5 6.249 E - 6 7.688 E - 7 9.604 E - 8 1.200 E - 8From Table3, we can numerically see that the h3-Richardson extrapolation algorithms are effective in obtaining high accuracy approximate solutions in which just the values at coarse grids and fine grids are applied. Furthermore, the a posteriori error estimate can be used to reckon the errors of the numerical solutions.


## 6. Conclusion

The mechanical quadrature methods provide a way to approximate the solutions of the nonlinear boundary value problems with high accuracy. The following conclusions can be drawn.These methods show that accuracy orders can be three and can be applied with the EAs to obtain higher accuracy. Furthermore, the larger scale of the problems, the more precision of the methods.To obtain the entry of discrete matrixes is very simple and straightforward, without any singular integrals.The problems, about Neumann boundary condition or nonlinear boundary condition in polygon elasticity, will be researched in our following work.



---

*Source: 6932164-2018-05-31.xml*