基于小波基函数的奇异积分求积算法*

杨焕枫,隆广庆

(南宁师范大学 数学与统计学院,广西 南宁 530100)

1 引言

工程和科学中的许多问题中都存在奇异积分,特别是在初值和边值问题中.例如在电磁波和声波散射问题中经常出现对数奇异核([1,2]);另一个著名的例子是二维亥姆霍兹方程,它可以通过格林公式直接转化为一个具有对数奇异核的积分方程([3]).奇异积分还出现在弱奇异积分算子的特征值问题中,比如文献[4]提出了一种基于矩阵压缩策略的弱奇异积分算子特征值问题的快速多尺度小波配置方法,最终得到了一个含奇异项的广义矩阵特征值问题;文献[5,6]则应用边界元的方法将微分方程转化为边界积分方程,得到了一个带有特殊结构的奇异核.因此,对此类积分的研究对工程和科学有着重要的意义.

在本文中,我们对多尺度逼近的背景下一类具有特殊结构的奇异核函数进行分析,建立适合使用多尺度小波基函数的快速求积方法.该方法能在一定程度上减少奇异积分的计算成本.我们定义这一类二重对数奇异积分为

其中[a,b]×[c,d]⊂[0,1]×[0,1],并且wij是小波基函数.

众所周知,在很多情况下奇异积分是无法准确计算的,因此产生了很多特殊方法用于处理奇异积分问题.例如,文献[7]提出了弱奇异积分的高斯型求积法则,并将该求积方法应用到第二类弱奇异Fredholm积分方程中,得到离散的乘积积分方法.文献[8]则主要是建立了一类在积分区域内或附近具有奇点的函数的数值积分的变阶复合求积公式,该数值方法具有指数收敛性.在文献[7,8]的基础上,文献[9]分别给出了具有多项式阶和指数阶的数值积分方法,并提出了误差控制策略,最后证明所给方法保证了收敛阶，减小了计算复杂度.之后,为了让求积法则适用于一般的奇异核类,文献[10]在文献[9]的研究基础上,提出适用于具有两个奇异点0和1的奇异函数的剖分方案的求积法则.对于弱奇异二重积分的计算,文献[11]给出了一个简单但有效的求积法则,文献[12,13]也对该方法进行了研究.此外,文献[14]还提出了一种计算具有奇异点{s-t=0}的二重奇异积分的自适应数值积分方法.由于该方法具有快速求解的优越性,在本文我们在该方法的基础上,给出一种求解具有奇异点{s-t=0,±1}的二重奇异积分的自适应数值积分方法.

近年来,乘积积分方法已广泛应用于计算奇异积分.例如,在文献[15,16]中用于离散奇异积分算子.而文献[5,6]为奇异核开发了一种特殊的乘积积分方法,利用其核函数的特殊结构能够开发一种比文献[10]中使用的高斯求积法则更有效的方法.为此,他们首先采用了一种技术,该技术在文献[17,18]中用于建立求解具有非线性边界条件的拉普拉斯方程和用于建立求解修正的亥姆霍兹方程的Fourier-Galerkin方法,对于矩阵中涉及的弱奇异积分的计算,他们利用乘积积分方法,即将弱奇异核分解为两个核的和,其中一个是携带主奇点,具有简单结构特征的弱奇异核,另一个是光滑核.简单奇异项能够在多尺度基下通过显式公式准确计算,而光滑项通过高精度求积方法计算,例如复合高斯勒让德求积法则.

在本文的第二节中,我们主要介绍以往计算奇异积分的数值积分方法.在第三节,我们重点对小波基函数进行数值处理,并给出线性小波基函数的构造方法.在第四节,针对一类对数奇异核,我们首先给出一种自适应数值积分方法来计算对数二重奇异积分,然后再提出一种能够准确计算出前面的二重对数奇异积分的乘积积分方法.最后一节给出的数值结果验证了本文提出的方法的有效性.

2 奇异积分的几种数值积分方法

在本节中,我们主要回顾几个计算弱奇异积分的求积法则,为后面验证算法的有效性做准备.

2.1 高斯勒让德型数值积分方法

我们先回顾文献[7]提出的弱奇异积分的高斯型方法.记Zm∶={0,1,…,m-1}.对于一个固定的正整数k,假设h∈C2k(0,1]且存在一个正常数c使得

|h(2k)(t)|≤c-σ-2k,t∈(0,1],σ∈[0,1).

现在考虑积分

使得Ij∶=[tj,tj+1],j∈Zm，是I∶=[0,1]的一部分.令

现在可以利用

2.2 多项式阶数值积分方法

假设h(s,t)=K(s,t)wij(t),其中核函数K(s,t)在s处具有奇异性,wij是小波基函数.同时还假设h具有下面的性质.

(I)Sij∶=supp(h)是I的一个子区间;

(II) 存在有限个点集π(h)∶={sj∶j-1∈Zm'-1}使得h∈C2k(I({s}∪π(h)));

(III) 存在一个正常数θ'使得|h(2k)(t)|≤θ'|t-s|-(σ+2k),t∈I({s}∪π(h)).

为了计算积分

选取与奇异点s相关联的两个节点集合

中的元素,并记它们为q'=q0

Π(h)∶={Qα∶=[qα,qα+1]∶α∈Zm'}.

为插值点的k-1次拉格朗日插值多项式.类似地我们可以利用I(Sk)去逼近I(h),再根据文献[9]中引理3.1，得到|Ek(f)|=O(m-2k).

可以看到,上面积分区间的剖分方案是适用于具有一个奇异点0的奇异函数.因此,为了构造一个更加适合奇异函数的数值积分,文献[10]在上面方法的基础上设计了适用于具有两个奇异点0和1的奇异函数的剖分方案.由于其剖分点为

因此,我们可以将该方案用到上述方法中去.

2.3 指数阶数值积分方法

同样考虑计算积分I(h),为此,选择一个新的区间剖分方案,即对于任意的γ∈(0,1),选取m+1个点

t0=0,tj=γm-j,j=1,2,…,m.

同样,可以设计适用于具有两个奇异点0和1的奇异函数的剖分方案.令

t0=0,tj=γm-j,j=1,2,…,m,以及tj=1-t2m-j,j=m+1,m+2,…,2m.

自然地也可以将这些剖分点应用到上述方法中去.

2.4 弱奇异二重积分的数值积分方法

接下来,简单回顾文献[11]中弱奇异二重积分的求积方法.考虑积分

其中D∶=[a,b]×[c,d],函数f,g都是光滑函数.

为了有效处理核函数在{(s,t):s-t=0,±1}上的奇异性,进行变量替换,令

则积分I等价于

并且

α(ξ)∶=max{2a-ξ,2c+ξ}和β(ξ)∶=min{2b-ξ,2d+ξ},

此时F(ξ,η)在点ξ=0,±1上具有奇异性,积分区域Ω∶={(ξ,η)∶|ξ|+|η-1|≤1}.

给定一个m>0,令

使用下面的4m+1个剖分点将区间[-1,1]划分为4m个子区间,即

-1+xj,-xj,xj,1-xj,j=0,1,…,m.

我们按照递增的顺序对这些剖分点重新排列,并记为ξi,i=0,1,…,4m,则

其中

需要注意的是,这里端点为奇异点0或±1的子区间被排除在计算之外.

在积分域为D=[a,b]×[c,d]的一般情况下,使用上述位于区间[a-d,b-c]内的剖分点,并重复上述过程.在这种情况下,还将端点a-d和b-c以及Ω的其他两个顶点的ξ坐标添加到剖分点中.

其中Δξ=(β(ξ)-α(ξ))/mξ,而ui和wi则分别对应勒让德多项式在[-1,1]上的零点和零点所对应的权重.

3 小波基函数的数值处理

我们知道,小波基函数不仅具有分片光滑性,而且还拥有有局部紧支集的良好性质,这驱使我们发展一种基于小波基函数的弱奇异积分的乘积积分方法,从而在实际应用中能大大地减少计算奇异积分的时间.接下来我们将对小波基函数进行数值处理,为我们下一节建立乘积积分方法做准备.

对于任意的i=0,j∈Zw(i),w0j是一个r次多项式.因此,我们可以把w0j写成下面的形式：

对于任意的i>0,j∈Zw(i),小波基函数wij是分片r次多项式.根据基函数wij的构造,对于任意的(i,j)∈Un,支撑集Sij可以被剖分成μ个子区间,即

Ωκ=(aκ,bκ),κ∈Zμ,

且在每个子区间上wij是r次多项式.因此,我们可以把wij写成

下面我们给出小波基函数的一个具体例子,为了简单,本文只考虑线性小波基函数(参考文献[11]).假设μ=2且选择以j/2n,j=1,2,…,2n-1为节点的[0,1]上的分片线性多项式的空间Xn.我们选择空间X0的基函数为

以及W1的基函数为

而空间Wi(i≥2)的基函数可以通过下式递归生成:

4 奇异积分的快速求积算法

在本节中, 我们首先给出对数奇异积分的自适应数值积分方法, 最后针对目标核函数进行讨论分析, 建立一种高效且准确的乘积积分方法.

4.1 自适应数值积分方法

本小节主要考虑奇异积分

其中核函数K(s,t)在{(s,t)∶s-t=0,±1}上具有奇异性.我们给出的方法的主要思想是根据被积函数的奇性对积分区间设计一套剖分方案, 然后运用复合非均匀高斯勒让德积分公式计算被积函数在每个子区间上的积分， 选取的求积节点数目在每个子区间上不同.该方法具有指数阶收敛性.

为此, 我们首先需要作类似于2.4节中的变量替换, 在此我们不再赘述.接下来我们对区间[0,1]剖分：

t0=0,tj=γm-j,j=1,2,…,m, 以及tj=1-t2m-j,j=m+1,m+2,…,2m.

定义一个点集:

Πγ,m∶={-1+xj,-xj,xj,1-xj∶j∈Zm+1}，

令

T∶=(Πγ,m∩{a-d,b-c})∪{a-d,b-c,a-c,b-d}.

按照递增顺序重新排列集合T中的元素, 并记它们为p'=p0

Π∶={Pβ∶=[pβ,pβ+1]∶β∈Zm''}.

为插值点的kj-1次拉格朗日插值多项式.至此,我们可以使用I(Sk)去逼近I(h).而对于正常积分h(ξ), 我们同样采用在2.4节中所描述的方法.

4.2 乘积积分方法

在这一小节里, 我们重点解决基于小波基函数的二重奇异积分

其中,K(s,t)∶=log(|s-t‖s-t-1‖s-t+1|),Sij∶=supp(wij),Si'j'∶=supp(wi'j')，的求积问题.

根据在第3节对小波基函数的讨论我们定义以下特殊积分.对于γ∈Zr,α∈∑∶={-1,0,1}以及a,b∈[0,1]且a

引理4.1若γ∈Zr,α∈以及a,b∈[0,1]且a

证明作分部积分可得

对于右边第二项的积分有

使用引理4.1中的积分公式, 我们可以准确计算下面两个奇异积分：

和

在引理4.1的基础上, 我们可以发展一种能够准确求解基于小波基函数的二重奇异积分的求积方法.

对于小波基函数wi'j'(s), 我们同样可以写成下面的形式：

对于γ,η∈Zr,α∈Σ∶={-1,0,1}以及a,b,c,d∈[0,1]且a

引理4.2若γ,η∈Zr,α∈以及a,b,c,d∈[0,1]且a

证明由引理4.1以及二项式定理可知

再对上式进行积分, 我们即可证得结论.

使用引理4.2中的积分公式, 我们可以准确计算下面四个二重奇异积分：

5 数值实验

在这一节中,我们给出使用上一节的乘积积分方法的两个数值例子.我们首先根据乘积积分方法去求奇异积分的准确解,然后将该解用到传统的数值积分方法中去计算误差,并根据该误差计算出对应的收敛阶,以此来验证我们方法的有效性.同时我们使用k=2的高斯求积作为奇异积分的数值评估.由于所考虑的都是对数型奇点的核函数K(s,t),因此，我们为积分区间剖分方案选择参数q=5或γ=0.17.

对于数值积分方法,我们给出两个计算收敛阶的公式：

另外,还记计算奇异积分所需的秒数为CT，我们所有的实验都是在具有2.2GHz频率和4G运行内存的个人计算机中执行,并且使用Fortran语言进行编程.

例1考虑积分

其中s=0.7.

利用乘积积分方法的显示公式,我们可以算得该奇异积分的准确解

I=-0.509 072 604 871 587 324 597 313 488 200 213,

其执行时间<0.01秒.数值结果列在表1和表2中.

表1 多项式阶数值积分方法的数值结果

表2 指数阶数值积分方法的数值结果

从表1、表2可以清楚地看到, 乘积积分方法给出了一维对数弱奇异积分的准确解, 并且计算速度比传统的数值积分方法快得多.

例2考虑积分

类似于例1, 我们算得准确解I=-0.477 411 277 760 219 301 121 078 492 187 340(执行时间<0.01秒).数值结果列在表3和表4中.

表3 文献[11]中的数值积分方法的数值结果

表4 自适应数值积分方法的数值结果

从表3、表4可见，对于二维对数弱奇异积分的计算, 在同样的精度下, 本文所建立的自适应数值积分方法以及乘积积分方法都比文献[11]中的数值积分方法高效且容易在计算机上实现. 这说明了本文给出的自适应数值积分方法以及乘积积分方法是有效的.