# Extrapolation Algorithm for Computing Multiple Cauchy Principal Value Integrals

Authors: Jin Li; Yongling Cheng
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

In this paper, the computation of multiple (including two dimensional and three dimensional) Cauchy principal integral with generalized composite rectangle rule was discussed, with the density function approximated by the middle rectangle rule, while the singular kernel was analytically calculated. Based on the expansion of the density function, the asymptotic expansion formulae of error functional are obtained. A series is constructed to approach the singular point, then the extrapolation algorithm is presented, and the convergence rate is proved. At last, some numerical examples are presented to validate the theoretical analysis.


---

## Body

## 1. Introduction

With the development of boundary element methods, a lot of attentions had recently been paid to Cauchy principal value integrals (including one dimensional and two dimensional). In this paper, we pay our attention to certain two-dimensional Cauchy principal value integral of the form(1)If,t,s=∫L1∫L2fx,yx−ty−sdxdy,where ∫L1∫L2 denotes a Cauchy Principle value integral L1=a,b,L2=c,d, and t,s∈a,b×c,d the singular point.As an adequate tool in boundary element methods, Cauchy principal value integral equations [1–7] have been appeared in many engineering problems [8], such as fluid mechanics, acoustics, fracture mechanics, elasticity, and electromagnetic scattering problems for the modeling of many physical situations. Numerous works such as the Gaussian method [9], the transformation method [10], and the Newton-Cote methods [11–15] have been devoted to developing efficient quadrature formulae. In recent years, the superconvergence phenomenon of the hypersingular integral defined on interval and circle is studied in [16–19] and the superconvergence phenomenon of the Cauchy principal value integral is also investigated in [20–22].In [1], based on the Gauss quadrature of one-dimensional Cauchy principal value integral with weight function ωx=1−xαl+xβ,α,β>−1, generalized quadrature Gauss rule for two-dimensional Cauchy principal value integrals is presented. In [23], for the problem of computing 2D Cauchy principal value integrals of the form ∫SFP0,PdP, authors have constructed quadrature rules of the Gaussian type. In [3], for the two-dimensional Cauchy principal value integrals with respect to generalized smooth Jacobi weight functions, authors have considered product rules of the Gauss type for the numerical approximation of two-dimensional Cauchy principal value integrals. In [24], a principal value definition of the basic hypersingular integral in the fundamental integral equation for two-dimensional cracks in three-dimensional isotropic elasticity is proposed. In [25], the classical composite rectangle (midpoint) rule for the computation of two-dimensional singular integrals is discussed, with the error functional of the rectangle rule for computing two-dimensional singular integrals, and the local coordinate of certain point and the convergence results Oh2 are obtained. In [26], the classical composite trapezoidal rule for the computation of two-dimensional singular integrals is presented and the convergence results Oh2 is the same as the Riemann integral convergence rate at a certain point of the classical composite trapezoidal rule. In [27], the Cauchy principal value integral is approximated by the zeros of the Chebyshev polynomials of the first kind. In [28], a twelve-point cubature formula has been constructed for the numerical evaluation of two-dimensional real Cauchy principal value integrals. In [29], Gaussian quadrature rules for the evaluation of Cauchy principal values of integrals are considered and their relation with Gauss-Legendre formulae is also studied. In [30], the density function fx,y was presented as(2)fx,y−ft,s=f1x,y+f2x,y+f3x,y,where(3)f1x,y=fx,s−ft,s,f2x,y=ft,y−ft,s,f3x,y=fx,y−fx,s−ft,y+ft,s.Then, (1) can be written as (4)∫L1∫L2fx,yx−ty−sdxdy=∫L1dyy−s∫L2f1x,yx−tdx+∫L1dxx−t∫L2f2x,yy−sdy+∫L1∫L2f3x,yx−ty−sdxdy+ft,s∫L1dxx−t∫L2dyy−s,which is the method with subtraction of singularity while the singular part was calculated analytically.In [31, 32], Lyness studied the Euler–Maclaurin expansion technique for the evaluation of two-dimensional Cauchy principal integrals. The extrapolation methods based on Hadamard finite-part integral definition with the trapezoidal rule for the computation of hypersingular integrals on interval and in circle in boundary element methods are presented in [33, 34], respectively. The extrapolation methods based on definition of subtraction of singularity are presented in [35, 36], respectively.In this paper, we focus on certain kinds of two-dimensional Cauchy principal value integrals which have not been extensively studied. It is the aim of this paper to investigate the error expansion of the rectangle rule and, in particular, to derive the error estimates. We examine the convergence property of the rectangle rule for certain kinds of the two-dimensional Cauchy principal value integrals and generalize the abovementioned one-dimensional convergence results to cover this new situation. Moreover, we give an error expansion of the corresponding remainder when the density functionfx,y belongs to C∞a,b×c,d. As the special function of the error expansion equals zero, we get the superconvergence phenomenon, i.e., when the singular point coincides with some a priori known point, the convergence rate is higher than what is globally possible. Then, the extrapolation algorithm is presented and the convergence rate of an extrapolation error and a posterior error is proved.The rest of this paper is organized as follows. In Section2, after introducing some basic formulae of the rectangle rule, we present the main results. In Section 3, we finished the proof. Finally, several numerical examples are provided to validate our analysis.


## 2. Main Result

Leta=x0<x1<⋯<xn−1<xn=b and c=y0<y1<⋯<ym−1<ym=d be a uniform partition of the area a,b×c,d with mesh size hx=b−a/n and hy=d−c/m. In order to simplify our analysis in the following, we set h=hx=hy=b−a/n=d−c/n; it is not difficult to extend our analysis to the case of quasi-uniform meshes.We setx^i=xi+h/2,y^j=yj+h/2,i,j=0,1,…,n−1, and define fCx,y as the piecewise constant interpolation for fx,y:(5)fCx,y=fx^i,y^j,i,j=0,1,…,n−1,and then we obtain(6)Inf,t,s=∫ab∫cdfCx,yx−ty−sdxdy=∑i=0n−1∑j=0n−1fx^i,y^j∫xixi+1∫yjyj+1dxdyx−ty−s=∑i=0n−1∑j=0n−1fx^i,y^jωi,jt,s=If,t,s−Enf,t,s,where Enf,t,s denotes the error functional and the Cotes coefficient is(7)ωi,jt,s=∫xixi+1∫yjyj+1dxdyx−ty−s=∫xixi+1dxx−t∫yjyj+1dyy−s=logxi+1−txi−tlogyj+1−syj−s.Define the following linear transformation:(8)x=x^iτ≔τ+1xi+1−xi2+xi,τ∈−1,1,y=y^jτ≔ξ+1yj+1−yj2+yj,ξ∈−1,1,from the reference element −1,1 to the subintervals xi,xi+1 and yj,yj+1.Let(9)ϕτt,s=−12∫−11∫−11ττ−tξ−sdτdξ,t<1,s<1,−12∫−11∫−11ττ−tξ−sdτdξ,t>1,s>1,(10)ϕξt,s=−12∫−11∫−11ξτ−tξ−sdτdξ,t<1,s<1,−12∫−11∫−11ξτ−tξ−sdτdξ,t>1,s>1,which can be expressed as the Legendre function of second kinds in [37]. We also set J≔x,y∈R2:x≠1andy≠1, and the operator W:CJ⟶C−1,1×−1,1 can be defined by(11)Wfτ,ξ≔fτ,ξ+∑i=1∞∑j=1∞f2i+τ,2j+ξ+∑i=1∞∑j=1∞f−2i+τ,−2j+ξ.Then, we have(12)S0xτ,ξ≔ϕxτ,ξ+∑i=1∞∑j=1∞ϕx2i+τ,2j+τ+∑i=1∞∑j=1∞ϕx−2i+τ,−2j+τ,(13)S0yτ,ξ≔ϕyτ,ξ+∑i=1∞∑j=1∞ϕy2i+τ,2j+τ+∑i=1∞∑j=1∞ϕy−2i+τ,−2j+τ.The superconvergence results of constant rectangle rules are given in the following.Theorem 1 (Theorem 2.2 in [20]).
LetS0xτ,ξandS0yτ,ξ be defined by (12) and (13), respectively. Assume that fx,y∈C2a,b×c,d and t≠xi,s≠yj, for any i,j=0,1,…,n. For the middle rectangle rule (6), there exists a positive constant C, independent of h, t,ands, such that(14)Enf,t,s=fyt,sS0xτ,ξh+fxt,sS0yτ,ξh+ℛnt,s,where t=xk+1+τ/2hands=yl+1+τ/2hk,l=0,1,…,n−1,τ∈−1,1.(15)ℛnt,s≤Clnh+lnγτ2+lnhηt+ηsh2,where(16)γτ=γξ=min0≤i≤n−1t−xih=min0≤j≤n−1s−yjh=1−τ2,(17)ηt=1t−ab−t,ηs=1s−cd−s.From Theorem1, when the special functions S0xτ,ξ=S0yτ,ξ=0, we get the superconvergence phenomenon. One can see that the superconvergence rate of the (composite) rectangle rule at certain points is one order higher than their global convergence rate which is the same as the Riemann integral.In the following theorem, we present the asymptotic expansion formulae of the error functional. The proof of this theorem will be given in the next part.Theorem 2.
Assume thatfx,y∈C∞a,b×c,d. For the rectangle rule Inf,t,s defined in 2.2, there exist certain constant ai,jτ,ξ, independent of h, t,ands, such that(18)Enf,t,s=∑l=1∞hl2l−1∑i,j=0,1,…,li+j=lfi,jt,sai,jτ,ξ,l≥1,where t=xk+1+τ/2hands=yl+1+ξ/2h,k,l=0,1,…,n−1,τ,ξ∈−1,1.


## 3. Proof of Theorem 2

In this section, we study the asymptotic expansion of the composite rectangle rule for the multiple (two dimensional) Cauchy principle value integrals.


### 3.1. Preliminaries

In the following analysis,C will denote a generic constant that is independent of h and s, and it may have different values in different places.Lemma 1.
Assumet,s∈xk,xk+1×yl,yl+1 and let ai=2t−xi/h−1,bj=2s−yj/h−1,0≤i,j≤n−1. Then, we have(19)ϕxai,bj=−1h∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy,i=k,j=l,−1h∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy,i≠k,j≠l,(20)ϕyai,bj=−1h∫xixi+1∫yjyj+1y−y^jx−ty−sdxdy,i=k,j=l,−1h∫xixi+1∫yjyj+1y−y^jx−ty−sdxdy,i≠k,j≠l.Proof.
In the case ofi=kandj=l, we have(21)∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy=∫xixi+1x−x^ix−tdx∫yjyj+11y−sdy=limε⟶0∫xit−ε+∫t+εxi+1x−x^ix−tdxlimε⟶0∫yjs−ε+∫s+εyj+11y−sdy=h2limε⟶0∫−1ai−2ε/h+∫ai+2ε/h1ττ−aidτlimε⟶0∫−1bj−2ε/h+∫bj+2ε/h11ξ−bjdξ=h2∫−11∫−11ττ−aiξ−bjdτdξ=−hϕxai,bj.[/lcmath]
Then, we prove the first identity in (19). The second identity in (19) can be obtained by applying the approach to the correspondent Riemann integral. The two identities in (20) can be obtained similarly as (19).Proof of Theorem 2.
By Taylor expansion, takingfCx,y and the point x,y, we obtain(22)∑i=0,i≠kn−1∑j=0,j≠ln−1∫xixi+1∫yjyj+1fx,y−fCx,yx−ty−sdxdy=−∑i=0,i≠kn−1∑j=0,j≠ln−1∫xixi+1∫yjyj+1x^i−xfxx,y+y^j−yfyx,yx−ty−sdxdy−∑i=0,i≠kn−1∑j=0,j≠ln−1∫xixi+1∫yjyj+1x^i−x2fxxx,y+y^j−y2fyyx,y+2x^i−xy^j−yfxyx,y2!x−ty−sdxdy−∑i=0,i≠kn−1∑j=0,j≠ln−1∫xixi+1∫yjyj+1∑p=03C3px^i−xpy^j−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∑i=0,i≠kn−1∑j=0,j≠ln−1∫xixi+1∫yjyj+1∑p=0mCmpx^i−xpy^j−ym−p∂m/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯.
For the casei=kandj=l, we have(23)∫xkxk+1∫ylyl+1fx,y−fCx,yx−ty−sdxdy=−∫xkxk+1∫ylyl+1x^k−xfxx,y+y^l−yfyx,yx−ty−sdxdy−∫xkxk+1∫ylyl+1x^k−x2fxxx,y+y^l−y2fyyx,y+2x^k−xy^l−yfxyx,y2!x−ty−sdxdy−∫xkxk+1∫ylyl+1∑p=03C3px^k−xpy^l−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∫xkxk+1∫ylyl+1∑p=0mCmpx^k−xpy^l−ym−p∂m/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯.
Putting (22) and (23) together, we have(24)∫ab∫cdfx,yx−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1fCx,yx−ty−sdxdy=∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1fx,yx−ty−sdxdy−∑i=0n−1∑j=0n−1∫xkxk+1∫ylyl+1fCx,yx−ty−sdxdy=∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1fx,y−fCx,yx−ty−sdxdy=−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−xfxx,y+y^j−yfyx,yx−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−x2fxxx,y+y^j−y2fyyx,y+2x^i−xy^j−yfxyx,y2!x−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=03C3px^i−xpy^j−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=0mCmpx^i−xpy^j−ym−p∂m/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯=fxt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y−y^jx−ty−sdxdy+fyt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−x2fxxt,s+y^j−y2fyyt,s+2x^i−xy^j−yfxyt,s2!x−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=03C3px^i−xpy^j−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=0mCmpx^i−xpy^j−ym−p∂m/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxt,sl!y−y^jx−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyt,sl!x−x^ix−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxxt,sl!x^i−x22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyyt,sl!y^j−y22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxyt,sl!2x^i−xy^j−y2!x−ty−sdxdy=fxt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y−y^jx−ty−sdxdy+fyt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy−fxxt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−x22!x−ty−sdxdy−fyyt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y^j−y22!x−ty−sdxdy−2fxyt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−xy^j−y2!x−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=03C3px^i−xpy^j−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=0mCmpx^i−xpy^j−ym−p∂3/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxt,sl!y−y^jx−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyt,sl!x−x^ix−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxxt,sl!x^i−x22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyyt,sl!y^j−y22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxyt,sl!2x^i−xy^j−y2!x−ty−sdxdy=fxt,shS0xτ,ξ+fyt,shS0yτ,ξ−12fxxt,sh2Sxxτ,ξ−12fyyt,sh2Syyτ,ξ−fxyt,sh2Syyτ,ξ−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=03C3px^i−xpy^j−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=0mCmpx^i−xpy^j−ym−p∂m/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxt,sl!y−y^jx−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyt,sl!x−x^ix−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxxt,sl!x^i−x22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyyt,sl!y^j−y22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxyt,sl!2x^i−xy^j−y2!x−ty−sdxdy=∑l=1∞hl2l−1∑i,j=0,1,…,li+j=lfi,jt,sai,jτ,ξ,l≥1,where we have used the Taylor expansion of fx,y at the point t,s:(25)fxx,y=fxt,s+∑l=1∞x−t∂/∂x+y−s∂/∂xlfxt,sl!,fyx,y=fyt,s+∑l=1∞x−t∂/∂t+y−s∂/∂slfyt,sl!.
Here, we have(26)S0xτ,ξ=1hlimn⟶∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y−y^jx−ty−sdxdy=1hlimn⟶∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1ττ−aiξ−bjdτdξ=limn⟶∞ϕxai,bj+∑i=1m∑j=1mϕx2j+τ,2j+ξ+∑i=0n−m∑j=1n−mϕx−2j+τ,−2j+ξ=ϕxτ,ξ+∑i=1∞∑j=1∞ϕx2i+τ,2j+ξ+∑i=1∞∑j=1∞ϕx−2i+τ,−2j+ξ.
Similarly, we can have(27)S0yτ,ξ=1hlimn⟶∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy=ϕyτ,ξ+∑i=1∞∑j=1∞ϕy2i+τ,2j+ξ+∑i=1∞∑j=1∞ϕy−2i+τ,−2j+ξ,(28)Sxxτ,ξ=1h2∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−x2x−ty−sdxdy=ϕxxτ,ξ+∑i=1∞∑j=1∞ϕxx2i+τ,2j+ξ+∑i=1∞∑j=1∞ϕxx−2i+τ,−2j+ξ,(29)Syyτ,ξ=1h2∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y^j−y2x−ty−sdxdy=ϕyyτ,ξ+∑i=1∞∑j=1∞ϕyy2i+τ,2j+ξ+∑i=1∞∑j=1∞ϕyy−2i+τ,−2j+ξ,(30)Sxyτ,ξ=1h2∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−xy^j−yx−ty−sdxdy=ϕxyτ,ξ+∑i=1∞∑j=1∞ϕxy2i+τ,2j+ξ+∑i=1∞∑j=1∞ϕxy−2i+τ,−2j+ξ,where(31)ϕxxai,bj=−1h2∫xixi+1∫yjyj+1x−x^i2x−ty−sdxdy,i=k,j=l,−1h2∫xixi+1∫yjyj+1x−x^i2x−ty−sdxdy,i≠k,j≠l,(32)ϕyyai,bj=−1h2∫xixi+1∫yjyj+1y−y^j2x−ty−sdxdy,i=k,j=l,−1h2∫xixi+1∫yjyj+1y−y^j2x−ty−sdxdy,i≠k,j≠l,(33)ϕxyai,bj=−1h2∫xixi+1∫yjyj+1x−x^iy−y^jx−ty−sdxdy,i=k,j=l,−1h2∫xixi+1∫yjyj+1x−x^iy−y^jx−ty−sdxdy,i≠k,j≠l.
As for the part of(34)∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y−y^jkx−x^ildxdy,which can be considered as the error estimate of the rectangle rule for the definite integral ∫ab∫cdx−tiy−sjdxdy,i,j≥1. Obviously, by the theorem ([38], Theorem D.2.6 P470), it can be expanded by the Euler–Maclaurin expansions, and we have(35)Enif,h=∫ab∫cdx−tiy−sjdxdy+∑k=1∞c2kh2k.
It is easy to see that there is no relation with the singular pointt,s and the Proof of Theorem 2 is completed.
Here, we can easily get thata01τ,ξ=S0yτ,ξanda10τ,ξ=S0xτ,ξ and a0,2τ,ξ=Syyτ,ξ,a2,0τ,ξ=Sxxτ,ξ,anda1,1τ,ξ=Sxyτ,ξ+Syxτ,ξ. For the one-dimensional Cauchy principal integral, we get the superconvergence point as ±2/3. As we have proved above, when a01τ,ξ=S0yτ,ξ=a10τ,ξ=S0xτ,ξ=0, then we get the superconvergence point of certain two-dimensional ±2/3,±2/3 which is the same as reference [20].
It is not difficult to extend our analysis to the case of multidimensional Cauchy principal value integrals. Letft1,t2,…,tk∈C∞ai,bikai=xi0<xi1<⋯<xi,n−1<xin=bi be a uniform partition of the area ai,bik with mesh size hkx=bki−aki/n. We set x^ki=xki+h/2,i=0,1,…,n−1, and define fCxk1,xk2,…,xkn as the piecewise constant interpolation for fx,y:(36)fCxk1,xk2,…,xkn=fx^k1,…,x^kn.
Then, we obtain(37)Inf;t1,t2,…,tk=∫a1b1⋯∫akbkfCx1,x2,…,xkx1−t1⋯xk−tkdx1⋯dxk=∑i1=0n−1⋯∑ik=0n−1fx^k1,…,x^knωi1,…,ikt1,t2,…,tk=If;t1,t2,…,tk−Enf;t1,t2,…,tk,where Enf;t1,t2,…,tk denotes the error functional and the Cotes coefficient is(38)ωi1,…,ikt1,t2,…,tk=logx1,i+1−t1x1i−t1⋯logxk,i+1−tkxki−tk.Theorem 3.
Assume thatft1,t2,…,tk∈C∞a,bk. For the rectangle rule Inf;t1,t2,…,tk defined in (36), there exist certain constant ai1,…,ikτ1,…,τk, independent of h and t1,t2,…,tk, such that(39)Enf,t1,t2,…,tk=∑l=1∞hl2l−1∑i1,…,ik=0,1,…,li1,…,ik=lfi1,…,ikt1,t2,…,tkai1,…,ikτ1,…,τk,where l≥1,t1=x1k+1+τ1/2h,…,tk=xkl+1+τk/2h,i1,…,ik=0,1,…,n−1.


## 3.1. Preliminaries

In the following analysis,C will denote a generic constant that is independent of h and s, and it may have different values in different places.Lemma 1.
Assumet,s∈xk,xk+1×yl,yl+1 and let ai=2t−xi/h−1,bj=2s−yj/h−1,0≤i,j≤n−1. Then, we have(19)ϕxai,bj=−1h∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy,i=k,j=l,−1h∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy,i≠k,j≠l,(20)ϕyai,bj=−1h∫xixi+1∫yjyj+1y−y^jx−ty−sdxdy,i=k,j=l,−1h∫xixi+1∫yjyj+1y−y^jx−ty−sdxdy,i≠k,j≠l.Proof.
In the case ofi=kandj=l, we have(21)∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy=∫xixi+1x−x^ix−tdx∫yjyj+11y−sdy=limε⟶0∫xit−ε+∫t+εxi+1x−x^ix−tdxlimε⟶0∫yjs−ε+∫s+εyj+11y−sdy=h2limε⟶0∫−1ai−2ε/h+∫ai+2ε/h1ττ−aidτlimε⟶0∫−1bj−2ε/h+∫bj+2ε/h11ξ−bjdξ=h2∫−11∫−11ττ−aiξ−bjdτdξ=−hϕxai,bj.[/lcmath]
Then, we prove the first identity in (19). The second identity in (19) can be obtained by applying the approach to the correspondent Riemann integral. The two identities in (20) can be obtained similarly as (19).Proof of Theorem 2.
By Taylor expansion, takingfCx,y and the point x,y, we obtain(22)∑i=0,i≠kn−1∑j=0,j≠ln−1∫xixi+1∫yjyj+1fx,y−fCx,yx−ty−sdxdy=−∑i=0,i≠kn−1∑j=0,j≠ln−1∫xixi+1∫yjyj+1x^i−xfxx,y+y^j−yfyx,yx−ty−sdxdy−∑i=0,i≠kn−1∑j=0,j≠ln−1∫xixi+1∫yjyj+1x^i−x2fxxx,y+y^j−y2fyyx,y+2x^i−xy^j−yfxyx,y2!x−ty−sdxdy−∑i=0,i≠kn−1∑j=0,j≠ln−1∫xixi+1∫yjyj+1∑p=03C3px^i−xpy^j−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∑i=0,i≠kn−1∑j=0,j≠ln−1∫xixi+1∫yjyj+1∑p=0mCmpx^i−xpy^j−ym−p∂m/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯.
For the casei=kandj=l, we have(23)∫xkxk+1∫ylyl+1fx,y−fCx,yx−ty−sdxdy=−∫xkxk+1∫ylyl+1x^k−xfxx,y+y^l−yfyx,yx−ty−sdxdy−∫xkxk+1∫ylyl+1x^k−x2fxxx,y+y^l−y2fyyx,y+2x^k−xy^l−yfxyx,y2!x−ty−sdxdy−∫xkxk+1∫ylyl+1∑p=03C3px^k−xpy^l−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∫xkxk+1∫ylyl+1∑p=0mCmpx^k−xpy^l−ym−p∂m/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯.
Putting (22) and (23) together, we have(24)∫ab∫cdfx,yx−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1fCx,yx−ty−sdxdy=∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1fx,yx−ty−sdxdy−∑i=0n−1∑j=0n−1∫xkxk+1∫ylyl+1fCx,yx−ty−sdxdy=∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1fx,y−fCx,yx−ty−sdxdy=−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−xfxx,y+y^j−yfyx,yx−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−x2fxxx,y+y^j−y2fyyx,y+2x^i−xy^j−yfxyx,y2!x−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=03C3px^i−xpy^j−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=0mCmpx^i−xpy^j−ym−p∂m/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯=fxt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y−y^jx−ty−sdxdy+fyt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−x2fxxt,s+y^j−y2fyyt,s+2x^i−xy^j−yfxyt,s2!x−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=03C3px^i−xpy^j−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=0mCmpx^i−xpy^j−ym−p∂m/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxt,sl!y−y^jx−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyt,sl!x−x^ix−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxxt,sl!x^i−x22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyyt,sl!y^j−y22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxyt,sl!2x^i−xy^j−y2!x−ty−sdxdy=fxt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y−y^jx−ty−sdxdy+fyt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy−fxxt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−x22!x−ty−sdxdy−fyyt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y^j−y22!x−ty−sdxdy−2fxyt,s∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−xy^j−y2!x−ty−sdxdy−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=03C3px^i−xpy^j−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=0mCmpx^i−xpy^j−ym−p∂3/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxt,sl!y−y^jx−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyt,sl!x−x^ix−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxxt,sl!x^i−x22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyyt,sl!y^j−y22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxyt,sl!2x^i−xy^j−y2!x−ty−sdxdy=fxt,shS0xτ,ξ+fyt,shS0yτ,ξ−12fxxt,sh2Sxxτ,ξ−12fyyt,sh2Syyτ,ξ−fxyt,sh2Syyτ,ξ−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=03C3px^i−xpy^j−y3−p∂3/∂xp∂y3−pfx,y3!x−ty−sdxdy−⋯−∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1∑p=0mCmpx^i−xpy^j−ym−p∂m/∂xp∂ym−pfx,ym!x−ty−sdxdy−⋯−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxt,sl!y−y^jx−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyt,sl!x−x^ix−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxxt,sl!x^i−x22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfyyt,sl!y^j−y22!x−ty−sdxdy−∑l=1∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−t∂/∂x+y−s∂/∂xlfxyt,sl!2x^i−xy^j−y2!x−ty−sdxdy=∑l=1∞hl2l−1∑i,j=0,1,…,li+j=lfi,jt,sai,jτ,ξ,l≥1,where we have used the Taylor expansion of fx,y at the point t,s:(25)fxx,y=fxt,s+∑l=1∞x−t∂/∂x+y−s∂/∂xlfxt,sl!,fyx,y=fyt,s+∑l=1∞x−t∂/∂t+y−s∂/∂slfyt,sl!.
Here, we have(26)S0xτ,ξ=1hlimn⟶∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y−y^jx−ty−sdxdy=1hlimn⟶∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1ττ−aiξ−bjdτdξ=limn⟶∞ϕxai,bj+∑i=1m∑j=1mϕx2j+τ,2j+ξ+∑i=0n−m∑j=1n−mϕx−2j+τ,−2j+ξ=ϕxτ,ξ+∑i=1∞∑j=1∞ϕx2i+τ,2j+ξ+∑i=1∞∑j=1∞ϕx−2i+τ,−2j+ξ.
Similarly, we can have(27)S0yτ,ξ=1hlimn⟶∞∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x−x^ix−ty−sdxdy=ϕyτ,ξ+∑i=1∞∑j=1∞ϕy2i+τ,2j+ξ+∑i=1∞∑j=1∞ϕy−2i+τ,−2j+ξ,(28)Sxxτ,ξ=1h2∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−x2x−ty−sdxdy=ϕxxτ,ξ+∑i=1∞∑j=1∞ϕxx2i+τ,2j+ξ+∑i=1∞∑j=1∞ϕxx−2i+τ,−2j+ξ,(29)Syyτ,ξ=1h2∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y^j−y2x−ty−sdxdy=ϕyyτ,ξ+∑i=1∞∑j=1∞ϕyy2i+τ,2j+ξ+∑i=1∞∑j=1∞ϕyy−2i+τ,−2j+ξ,(30)Sxyτ,ξ=1h2∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1x^i−xy^j−yx−ty−sdxdy=ϕxyτ,ξ+∑i=1∞∑j=1∞ϕxy2i+τ,2j+ξ+∑i=1∞∑j=1∞ϕxy−2i+τ,−2j+ξ,where(31)ϕxxai,bj=−1h2∫xixi+1∫yjyj+1x−x^i2x−ty−sdxdy,i=k,j=l,−1h2∫xixi+1∫yjyj+1x−x^i2x−ty−sdxdy,i≠k,j≠l,(32)ϕyyai,bj=−1h2∫xixi+1∫yjyj+1y−y^j2x−ty−sdxdy,i=k,j=l,−1h2∫xixi+1∫yjyj+1y−y^j2x−ty−sdxdy,i≠k,j≠l,(33)ϕxyai,bj=−1h2∫xixi+1∫yjyj+1x−x^iy−y^jx−ty−sdxdy,i=k,j=l,−1h2∫xixi+1∫yjyj+1x−x^iy−y^jx−ty−sdxdy,i≠k,j≠l.
As for the part of(34)∑i=0n−1∑j=0n−1∫xixi+1∫yjyj+1y−y^jkx−x^ildxdy,which can be considered as the error estimate of the rectangle rule for the definite integral ∫ab∫cdx−tiy−sjdxdy,i,j≥1. Obviously, by the theorem ([38], Theorem D.2.6 P470), it can be expanded by the Euler–Maclaurin expansions, and we have(35)Enif,h=∫ab∫cdx−tiy−sjdxdy+∑k=1∞c2kh2k.
It is easy to see that there is no relation with the singular pointt,s and the Proof of Theorem 2 is completed.
Here, we can easily get thata01τ,ξ=S0yτ,ξanda10τ,ξ=S0xτ,ξ and a0,2τ,ξ=Syyτ,ξ,a2,0τ,ξ=Sxxτ,ξ,anda1,1τ,ξ=Sxyτ,ξ+Syxτ,ξ. For the one-dimensional Cauchy principal integral, we get the superconvergence point as ±2/3. As we have proved above, when a01τ,ξ=S0yτ,ξ=a10τ,ξ=S0xτ,ξ=0, then we get the superconvergence point of certain two-dimensional ±2/3,±2/3 which is the same as reference [20].
It is not difficult to extend our analysis to the case of multidimensional Cauchy principal value integrals. Letft1,t2,…,tk∈C∞ai,bikai=xi0<xi1<⋯<xi,n−1<xin=bi be a uniform partition of the area ai,bik with mesh size hkx=bki−aki/n. We set x^ki=xki+h/2,i=0,1,…,n−1, and define fCxk1,xk2,…,xkn as the piecewise constant interpolation for fx,y:(36)fCxk1,xk2,…,xkn=fx^k1,…,x^kn.
Then, we obtain(37)Inf;t1,t2,…,tk=∫a1b1⋯∫akbkfCx1,x2,…,xkx1−t1⋯xk−tkdx1⋯dxk=∑i1=0n−1⋯∑ik=0n−1fx^k1,…,x^knωi1,…,ikt1,t2,…,tk=If;t1,t2,…,tk−Enf;t1,t2,…,tk,where Enf;t1,t2,…,tk denotes the error functional and the Cotes coefficient is(38)ωi1,…,ikt1,t2,…,tk=logx1,i+1−t1x1i−t1⋯logxk,i+1−tkxki−tk.Theorem 3.
Assume thatft1,t2,…,tk∈C∞a,bk. For the rectangle rule Inf;t1,t2,…,tk defined in (36), there exist certain constant ai1,…,ikτ1,…,τk, independent of h and t1,t2,…,tk, such that(39)Enf,t1,t2,…,tk=∑l=1∞hl2l−1∑i1,…,ik=0,1,…,li1,…,ik=lfi1,…,ikt1,t2,…,tkai1,…,ikτ1,…,τk,where l≥1,t1=x1k+1+τ1/2h,…,tk=xkl+1+τk/2h,i1,…,ik=0,1,…,n−1.


## 4. Extrapolation Method

In the abovementioned sections, we have proved that the error functional of the middle rectangle rule has the following asymptotic expansion:(40)Enf,t,s=∑l=1∞hl2l−1∑i,j=0,1,…,li+j=lfi,jt,sai,jτ,ξ,l≥1.It is easy to see that the error functional depends on the value ofai,jτ,ξ.Now we present an algorithm for the givent,s. Assume there exist positive integers n01andn02 such that(41)m10≔n10t−ab−a,m20≔n20s−cd−c,are positive number. In order to simplify the deduction process, we assume that n10=m20. Firstly, we partition a,b×c,d into n01 equal subinterval to get a mesh denoted by Π1 with mesh size h1=b−a/n10=d−c/n20. Then, we refine Π1 to get mesh Π2 with mesh size h2=h1/2. In this way, we get a series of meshes Πjj=1,2,… in which Πj is refined from Πj−1 with the mesh size denoted by hj. The extrapolation scheme is presented in Table 1.Table 1
Extrapolation scheme ofTki.
T
h
1
=
T
1
1
T
h
2
=
T
1
2
T
2
1
T
h
3
=
T
1
3
T
2
2
T
3
1
T
h
4
=
T
1
4
T
2
3
T
3
2
T
4
1
T
h
5
=
T
1
5
T
2
4
T
3
3
T
4
2
T
5
1
⋮
⋮
⋮
⋮
⋮For a givenτ,ξ∈−1,1, define(42)ti=t+τ+12hi,i=1,2,…,sj=s+ξ+12hj,j=1,2,…,(43)Thi=Thi,hj=I2i−1n10,2j−1n20f,ti,sj.We present the following extrapolation algorithm.Step one:(44)ComputeT1j=Thj,j=1,…,m.Step two:(45)ComputeTij=Ti−1j+1+Ti−1j+1−Ti−1j2i−1−1,i=2,…,m,j=1,…,m−i.Theorem 4.
Under the asymptotic expansion of Theorem2, for a given τ,ξ and the series of meshes defined by (42), we have(46)If,t,s−Tij≤Chi,and a posteriori asymptotic error estimate is given by(47)Ti−1j+1−Ti−1j2i−1−1≤Chi−1,i≥2.Proof.
By the asymptotic expansion of (40), for given τ, we have(48)If,t,s−Thi=If,t,s−If,ti,sj+If,ti,sj−Thi=If,t,s−If,ti,sj+∑l=1∞hl2l−1∑i,j=0,1,…,li+j=lfi,jt,sai,jτ,ξ.
By the definition of Cauchy principal integral and (42), for the first two parts of (48) by Taylor expansion for If,ti,sj at the point t,s,(49)If,ti,sj=If;t,s+∑i,j=0,1,1Ii,jf,t,sτ+12hi+ξ+12hi+∑i,j=0,1,22Ii,jf,t,sτ+12hi+ξ+12hi2+⋯+∑i,j=0,1,…,llIi,jf,t,sτ+12hi+ξ+12hil+⋯.
Similarly we also expandfi,jti,sj at point t,s, and then we have(50)fi,jti,sj=fi,jt,s+∑k,l=0,1,1fi+k,j+lt,sτ+12hi+ξ+12hi+∑k,l=0,1,22fi+k,j+lt,s2!τ+12hi+ξ+12hi2+⋯+∑k,l=0,1,…,mmfi+k,j+lt,sl−i−1!τ+12hi+ξ+12hil−i−1+⋯.
Putting (48)–(50) together, we have(51)If,t,s−Thi=∑i=1∞bit,s,τh1i,where(52)bit,s,τ,ξ=∑i,j=0,1,…,llfi,jt,s∑k=1iai,jτ,ξ2kτ+12+ξ+12i−k1i−k!−τ+1i2ii!∑i,j=0,1,…,llIi,jf,t,s.
Then,bit,s,τ,ξ is a constant for a given τ and ξ. By (51), we can obtain(53)If,t,s−Thi+1=∑i=1∞bit,s,τh2i+1.
By equations (51) and (53), with hj=2hj+1, we easily have(54)If,t,s=2Thi+1−Thi+∑i=2∞bit,s,τ12i−1−1hij=T2j+∑i=2∞bit,s,τ12i−1−1hij,which means(55)If,t,s−T2j=∑i=2∞bit,s,τ12i−1−1hji,(56)T2j=2Thi+1−Thi.
Following the extrapolation process above, we can get the accuracyOh3, and we continue to use extrapolation process again and get Oh4. In this way, we finish the proof by mathematical induction.


## 5. Numerical Example

In this section, computational results are reported to confirm our theoretical analysis.Example 1.
Consider the Cauchy principal value integral(57)∫−11∫−11x3y3x−ty−sdxdy,t,s∈−1,1×−1,1,with the exact solution(58)t3log1−t1+t+2t2+23s3log1−s1+s+2s2+23.
We adopt the uniform meshes to examine the convergence rate of the rectangle rule with the dynamic point witht=xn/4+1+τh/2ands=yn/4+ξ+1h/2. From Table 2, we know that when the local coordinate of the singular point is ±2/3,±2/3, the quadrature reaches the convergence rate of Oh2; as for the nonsupersingular point, the convergence rate is Oh which agree with our theoretical analysis. For the case of t=xn/4+τ+1h/2ands=b−τ+1h/2 and t=a+τ+1h/2ands=b−τ+1h/2, Tables 3 and 4 show that there is no superconvergence phenomenon for the midrectangle rue which coincide with our theoretical analysis.Table 2
Errors of the rectangle rule witht=xn/4+τ+1h/2 and s=yn/4+ξ+1h/2.
n
0,0
2
/
3
,
2
/
3
−
2
/
3
,
−
2
/
3
2
/
3
,
−
2
/
3
−
2
/
3
,
2
/
3
64
5.8370e − 02
5.1411e − 03
−1.4502e − 03
1.9813e − 03
1.9813e − 03
128
3.1292e − 02
1.4307e − 03
−3.6303e − 04
5.5081e − 04
5.5081e − 04
256
1.6182e − 02
3.7674e − 04
−9.0779e − 05
1.4508e − 04
1.4508e − 04
512
8.2257e − 03
9.6619e − 05
−2.2695e − 05
3.7224e − 05
3.7224e − 05
1024
4.1467e − 03
2.4462e − 05
−5.6736e − 06
9.4270e − 06
9.4270e − 06
h
α
0.954
1.9288
1.9994
1.9289
1.9289Table 3
Errors of the rectangle rule witht=xn/4+τ+1h/2 and s=b−ξ+1h/2
n
0,0
2
/
3
,
2
/
3
−
2
/
3
,
−
2
/
3
2
/
3
,
−
2
/
3
−
2
/
3
,
2
/
3
64
−7.7845e − 02
−1.8981e − 02
−6.3590e − 02
−6.9936e − 02
−1.4985e − 02
128
−5.0142e − 02
−9.6972e − 03
−3.3326e − 02
−3.5663e − 02
−8.1438e − 03
256
−3.0714e − 02
−4.8309e − 03
−1.7118e − 02
−1.7883e − 02
−4.2928e − 03
512
−1.8182e − 02
−2.3911e − 03
−8.6906e − 03
−8.9255e − 03
−2.2176e − 03
1024
−1.0503e − 02
−1.1842e − 03
−4.3828e − 03
−4.4522e − 03
−1.1308e − 03
h
α
0.7225
1.0006
0.9647
0.9934
0.9320Table 4
Errors of the rectangle rule witht=a+τ+1h/2 and s=b−ξ+1h/2
n
0,0
2
/
3
,
2
/
3
−
2
/
3
,
−
2
/
3
2
/
3
,
−
2
/
3
−
2
/
3
,
2
/
3
64
−2.3835e − 01
−2.3005e − 02
−2.3005e − 02
−2.3036e − 02
−5.7772e − 02
128
−1.5647e − 01
−8.5913e − 03
−8.5913e − 03
1.2073e − 02
−4.0355e − 02
256
−9.8100e − 02
−2.9871e − 03
−2.9871e − 03
1.7480e − 02
−2.6800e − 02
512
−5.9436e − 02
−9.8785e − 04
−9.8785e − 04
1.4112e − 02
−1.7061e − 02
1024
−3.5070e − 02
−3.1487e − 04
−3.1487e − 04
9.5823e − 03
−1.0489e − 02
h
α
0.6912
1.5478
1.5939
0.3164
0.6154Example 2.
In this example, we still consider the Cauchy principal value integral(59)∫−11∫−11x3y6x−ty−sdxdy,t,s∈−1,1×−1,1,with the exact solution(60)t3log1−t1+t+2t2+232s5+23s3+25s+s6log1−s1+s.
Tables5 and 6 show extrapolation errors of the rectangle rule with the singular points s=xn/4+τ+1h/2andt=yn/4+τ+1h/2; the convergence rate of the rectangle rule is Oh,Oh2,Oh3,…, while a posterior error of the rectangle rule extrapolation convergence rate of the rectangle rule is Oh,Oh2,Oh3,…, with τ=ξ=0, which agree with our Theorem 4.
Tables7 and 8 show extrapolation errors of the rectangle rule with the singular points s=xn/4+τ+1h/2andt=yn/4+τ+1h/2, which show that the convergence rate of the rectangle rule is Oh,Oh2,Oh3,…, while a posterior error of the rectangle rule extrapolation convergence rate of the rectangle rule is Oh,Oh2,Oh3,…, with τ=ξ=2/3, which agrees with our Theorem 4.Table 5
Extrapolation errors of the rectangle rule withs=xn/4+τ+1h/2,t=yn/4+τ+1h/2,andτ=ξ=0.
h
−
extra
h
2
−
extra
h
3
−
extra
h
4
−
extra
16
−1.1103e − 01
32
−6.1657e − 02
−1.2288e − 02
64
−3.2546e − 02
−3.4355e − 03
−4.8471e − 04
128
−1.6726e − 02
−9.0658e − 04
−6.3594e − 05
−3.4340e − 06
256
−8.4795e − 03
−2.3270e − 04
−8.0783e − 06
−1.4752e − 07
512
−4.2692e − 03
−5.8938e − 05
−1.0154e − 06
−6.4608e − 09
1024
−2.1420e − 03
−1.4830e − 05
−1.2720e − 07
−3.0565e − 10Table 6
A posterior error of the rectangle rule withs=xn/4+τ+1h/2,t=yn/4+τ+1h/2,andτ=ξ=0.
h
−
extra
h
2
−
extra
h
3
−
extra
h
4
−
extra
16
32
−4.9369e − 02
64
−2.9111e − 02
−2.9508e − 03
128
−1.5820e − 02
−8.4299e − 04
−6.0160e − 05
256
−8.2468e − 03
−2.2463e − 04
−7.9308e − 06
−2.1910e − 07
512
−4.2103e − 03
−5.7922e − 05
−1.0090e − 06
−9.4039e − 09
1024
−2.1272e − 03
−1.4703e − 05
−1.2689e − 07
−4.1034e − 10Table 7
Extrapolation errors of the rectangle rule withs=xn/4+τ+1h/2,t=yn/4+τ+1h/2,andτ=ξ=2/3.
h
−
extra
h
2
−
extra
h
3
−
extra
h
4
−
extra
16
−2.3988e − 01
32
−1.1270e − 01
1.4482e − 02
64
−5.1968e − 02
8.7623e − 03
6.8556e − 03
128
−2.4458e − 02
3.0524e − 03
1.1491e − 03
3.3393e − 04
256
−1.1785e − 02
8.8738e − 04
1.6570e − 04
2.5214e − 05
512
−5.7733e − 03
2.3852e − 04
2.2231e − 05
1.7349e − 06
1024
−2.8558e − 03
6.1789e − 05
2.8785e − 06
1.1385e − 07Table 8
A posterior error of the rectangle rule withs=xn/4+τ+1h/2,t=yn/4+τ+1h/2,andτ=ξ=2/3.
h
−
extra
h
2
−
extra
h
3
−
extra
h
4
−
extra
16
32
−1.2718e − 01
64
−6.0730e − 02
1.9067e − 03
128
−2.7510e − 02
1.9033e − 03
1.4049e − 04
256
−1.2673e − 02
7.2168e − 04
8.1521e − 04
2.0581e − 05
512
−6.0118e − 03
2.1629e − 04
2.0496e − 05
1.5653e − 06
1024
−2.9175e − 03
5.8910e − 05
2.7647e − 06
1.0807e − 07Example 3.
Consider the Cauchy principal value integral(61)∫−11∫−11∫−11x3y6z3x−ty−sz−udxdydz,t,s,u∈−1,1×−1,1×−1,1,with the exact solution(62)t3log1−t1+t+2t2+232s5+23s3+25s+s6log1−s1+su3log1−u1+u+2u2+23.
Tables9 and 10 show extrapolation errors of the rectangle rule with the singular points t=xn/4+τ1+1h/2,s=xn/4+τ2+1h/2,andu=xn/4+τ3+1h/2; the convergence rate of rectangle rule is Oh,Oh2,Oh3,…, for the general convergence rate, while a posterior error of the rectangle rule extrapolation convergence rate of the rectangle rule is Oh,Oh2,Oh3,…, with τ1=τ2=τ3=0, which agrees with our Theorem 4.
Tables11 and 12 show extrapolation errors of the rectangle rule with the singular points t=xn/4+τ1+1h/2,s=xn/4+τ2+1h/2,andu=xn/4+τ3+1h/2; the convergence rate of the rectangle rule is Oh,Oh2,Oh3,…, for the general convergence rate, while a posterior error of the rectangle rule extrapolation convergence rate of the rectangle rule is Oh,Oh2,Oh3,…, with τ1=τ2=τ3=2/3, which agrees with our Theorem 4.Table 9
Extrapolation errors of the rectangle rule witht=xn/4+τ1+1h/2,s=xn/4+τ2+1h/2,u=xn/4+τ3+1h/2,s,τ1=τ2=τ3=0.
h
−
extra
h
2
−
extra
h
3
−
extra
h
4
−
extra
8
−2.2044e − 01
16
−1.3967e − 01
−5.8905e − 02
32
−7.9385e − 02
−1.9096e − 02
−5.8267e − 03
64
−4.2405e − 02
−5.4249e − 03
−8.6774e − 04
−1.5932e − 04
128
−2.1924e − 02
−1.4433e − 03
−1.1610e − 04
−8.7182e − 06
256
−1.1148e − 02
−3.7200e − 04
−1.4901e − 05
−4.4454e − 07Table 10
A posterior error of the rectangle rule witht=xn/4+τ1+1h/2,s=xn/4+τ2+1h/2,u=xn/4+τ3+1h/2,s,τ1=τ2=τ3=0.
h
−
extra
h
2
−
extra
h
3
−
extra
h
4
−
extra
8
16
−8.0769e − 02
32
−6.0289e − 02
−1.3270e − 02
64
−3.6980e − 02
−4.5571e − 03
−7.0842e − 04
128
−2.0481e − 02
−1.3272e − 03
−1.0738e − 04
−1.0040e − 05
256
−1.0776e − 02
−3.5710e − 04
−1.4456e − 05
−5.5157e − 07Table 11
Extrapolation errors of the rectangle rule witht=xn/4+τ1+1h/2,s=xn/4+τ2+1h/2,u=xn/4+τ3+1h/2,s,τ1=τ2=τ3=2/3.
h
−
extra
h
2
−
extra
h
3
−
extra
h
4
−
extra
8
−2.2044e − 01
16
−1.3967e − 01
−5.8905e − 02
32
−7.9385e − 02
−1.9096e − 02
−5.8267e − 03
64
−4.2405e − 02
−5.4249e − 03
−8.6774e − 04
−1.5932e − 04
128
−2.1924e − 02
−1.4433e − 03
−1.1610e − 04
−8.7182e − 06
256
−1.1148e − 02
−3.7200e − 04
−1.4901e − 05
−4.4454e − 07Table 12
A posterior error of the rectangle rule witht=xn/4+τ1+1h/2,s=xn/4+τ2+1h/2,u=xn/4+τ3+1h/2,s,τ1=τ2=τ3=2/3.
h
−
extra
h
2
−
extra
h
3
−
extra
h
4
−
extra
8
16
−8.0769e − 02
32
−6.0289e − 02
−1.3270e − 02
64
−3.6980e − 02
−4.5571e − 03
−7.0842e − 04
128
−2.0481e − 02
−1.3272e − 03
−1.0738e − 04
−1.0040e − 05
256
−1.0776e − 02
−3.5710e − 04
−1.4456e − 05
−5.5157e − 07


## 6. Conclusion

In this paper, we have shown, both theoretically and numerically, that the main part of error functional of the rectangle rule has the asymptotic expansion (18). Numerical experiment has shown that the special function ai,jτ,ξ has a big influence on the convergence rate. Moreover, by the extrapolation method, we not only obtain a high order of accuracy but also derive a posteriori error estimate conveniently. In fact, it is not difficult to extend our results to arbitrary (multidimensional) Cauchy principal value integral and the extrapolation methods can be similarly obtained.



---

*Source: 6340231-2020-03-11.xml*