塔形網(wǎng)格上的有理插值

許雨靜,趙前進

(安徽理工大學(xué) 數(shù)學(xué)與大數(shù)據(jù)學(xué)院,安徽 淮南 232001)

0 引 言

有理插值是用于解決非線性逼近問題不可或缺的方法.連分式因具有計算簡便、運算量少等特性,成為構(gòu)造有理插值算法的重要工具.近年來,在矩形網(wǎng)格、三角網(wǎng)格上構(gòu)造混合有理插值用以逼近二元和三元連續(xù)函數(shù)的算法越來越多[1-7],帶有附加條件的切觸有理插值[8-9]也很常見.沈曉明等將重心插值與Thiele連分式插值結(jié)合,得到了二元Barycentric-Thiele混合有理插值[7].這種插值方法既有重心插值的優(yōu)點,又有連分式的優(yōu)點,如通過選擇合適的權(quán)值,可以避免極點的產(chǎn)生,且計算量小等.文獻[10]證實,在特定條件下,連分式插值還具有保水平漸近線的性質(zhì).

現(xiàn)有塔型網(wǎng)格上的Thiele型分叉連分式插值算法[11]要定義初始點的值,這給計算過程帶來不便.針對這一情況,出現(xiàn)了不同塔型格式的Newton-Thiele型混合有理插值.但數(shù)值實驗發(fā)現(xiàn),這種算法在x=0處附近的插值效果并不理想.為彌補這一缺陷,本文將構(gòu)造新的Thiele型分叉連分式算法,并用實例驗證其逼近效果要優(yōu)于Newton-Thiele算法.

給定實數(shù)點集X={x0,x1,…,xm}∈[a,b]?,以及f(x)∈C1[a,b]是定義在區(qū)間上的連續(xù)實函數(shù),其中,當(dāng)xi≠xj時,有f(xi)≠f(xj).

(1)

其中,bk=φ[x0,x1,…,xk],φ[x0,x1,…,xk]稱為f(x)的k階逆差商,其定義為:

φ[xp]=f(xp),

(2)

(3)

(4)

T(x)稱為函數(shù)f(x)的Thiele型插值函數(shù).

1 連分式有理插值算法

將Thiele型連分式應(yīng)用于新的塔型網(wǎng)格上,對Thiele型連分式插值進行類似張量積的處理,從而在塔形網(wǎng)格上構(gòu)造分叉連分式插值,并給出遞推算法,以證明該插值的有效性.

(5)

其中,i=0,1,…,n.則有:

(6)

當(dāng)j、k=0,1,…,2(n-i)時,有:

(7)

(8)

其中,

(9)

(10)

式(5)～(8)為連分式有理插值函數(shù).

下面證明其有效性.

φ[xi,yj,zk]=f(xi,yj,zk),

(11)

(12)

(13)

(14)

(15)

(16)

(17)

(18)

定理1若定義1中的所有混合偏逆差商均存在:

其中,i=0,1,…,n;j、k=0,1,…,2(n-i).則有理插值函數(shù)滿足:

Rn(xi,yj,zk)=f(xi,yj,zk).

(19)

因此,

同理可證:

當(dāng)i=0,1,…,n時,有:

故由式(5)得:

定理1得證.

2 特征定理

引理1對i=0,1,…,n,令

(20)

則

(21)

(22)

證明證明方法同文獻[11]的引理4.5,令n=n-1即可.

(23)

degyPn(x,y,z)=degyQn(x,y,z)=2n(n+1)(n+2)/3,

(24)

degzPn(x,y,z)=degzQn(x,y,z)=2n(n+1)(n+2)/3.

(25)

證明式(23)可由Thiele型連分式插值最高次公式[12]得到.

下面證明式(24)和式(25).

當(dāng)n=1時,

顯然成立.

假設(shè)當(dāng)n=k時成立,即

degyPk(x,y,z)=degyQk(x,y,z)=k(k+1)(k+2)/3,
degzPk(x,y,z)=degzQk(x,y,z)=k(k+1)(k+2)/3.

下面證明當(dāng)n=k+1時,結(jié)論也滿足:

為證明該定理,不妨記:

其中,τ=k,k+1,…

其中,τ=k-1,k,…

不難看出,

且

因此,可得:

同理可證:

degzPk+1(x,y,z)=degzQk+1(x,y,z)=2(k+1)(k+2)(k+3)/3.

假設(shè)成立,該定理得證.

3 誤差分析

下面給出連分式有理插值函數(shù)的誤差估計.

定理3假設(shè)函數(shù)f(x,y,z)是包含在區(qū)間Ω=[a,b]×[c,d]×[e,f]上的三元函數(shù),且混合偏逆差商均存在,則對?(x,y,z)∈G,存在(θi,δj,i,ζi)∈G,(αi,βi)∈[c,d]×[e,f],使得

(26)

其中,Rn(x,y,z)=Pn(x,y,z)/Qn(x,y,z).

當(dāng)i=0,1,…,n;j、k=0,1,…,2(n-i)時,有:

(27)

θi∈I[x0,…,xi],δj,i∈I[y0,…,yj],ξi∈I[z0,…,z2(n-i),z],αi∈I[x0,…,xi],βi∈I[y0,…,y2(n-i),y],τ∈I[x0,…,xi,x],且I[x0,…,xi]是包含x0,…,xi的最小區(qū)間.

證明令En(x,y,z)=Qn(x,y,z)f(x,y,z)-Pn(x,y,z),(x,y,z)∈G.

En(xi,yj,zk)=0.

因此,基于張量積形式的Newton差商記為En[x0,…,xi;y0,…,yi;z0,…,zk],且En[x0,…,xi;y0,…,yi;z0,…,zk]=0.

利用展開式,可得:

4 數(shù)值舉例

下面給出塔型網(wǎng)格Newton-Thiele型混合有理插值算法與連分式有理插值算法的數(shù)值例子,并對結(jié)果進行對比與分析.

不妨記RT(x,y,z)為連分式有理插值函數(shù),RN(x,y,z)為Newton-Thiele型混合有理插值函數(shù),eT(x,y,z)、eN(x,y,z)分別為RT(x,y,z)、RN(x,y,z)與被插值函數(shù)的絕對誤差,即

eT(x,y,z)=|f(x,y,z)-RT(x,y,z)|;eN(x,y,z)=|f(x,y,z)-RN(x,y,z)|.

例1設(shè)函數(shù)f(x,y,z)=(x2+y2+z2)e-(x+y+z),且(x,y,z)∈{(xi,yj,zk)|i=0,1;j、k=0,1,…,2(1-i)},不妨令n=2,則給定的插值點見表1和表2.

表1 x0=0 初始表

表2 x1=1初始表

當(dāng)i=2時,x2=2,y0=1/2,z0=1/3;f(x1,y0,z0)=0.256 5.

由(5)～(18)式,可得:

當(dāng)i=0,1,2時,有:

當(dāng)i=0,j、k=0,1,2,3,4時,有:

其中,

同理,當(dāng)i=1;j、k=0,1,2時,可求bi,j(z),ci,k(y).

當(dāng)i=2;j、k=0時,有:

通過計算,可得連分式有理插值系數(shù),見表3、表4和表5.

表4 RT(x,y,z)、RN(x,y,z)系數(shù)表

表5 RT(x,y,z)插值c1,j(y),b1,j(z)系數(shù)表

下面計算Newton-Thiele型混合有理插值系數(shù).

通過計算得,當(dāng)i=0時,插值系數(shù)見表3和表4.

當(dāng)i=1時,插值系數(shù)見表6.

表6 RN(x,y,z)插值c1,j(y),b1,j(z)系數(shù)表

當(dāng)i=2時,插值系數(shù)為:

通過大量的實驗對比發(fā)現(xiàn),當(dāng)x的取值趨于0,y,z取定義域內(nèi)的任意值時,連分式有理插值算法的誤差明顯小于Newton-Thiele型混合有理插值算法的誤差,且插值效果較好,表7僅代表一小部分?jǐn)?shù)據(jù).

圖1 x=0.01,y=0.8,z∈[1/3,1]函數(shù)圖像

圖2 x=0.001,y=0.8,z∈[1/3,1]函數(shù)圖像

5 結(jié) 論

Newton插值定義于差商上,而Thiele型連分式定義于逆差商上,這就會導(dǎo)致兩者在插值問題上出現(xiàn)不同.文中數(shù)值例子證明,在某些點的周圍,連分式有理插值要優(yōu)于Newton-Thiele型混合有理插值.但在整個定義域內(nèi),前者與后者的逼近效果不相上下.因此,后續(xù)研究的方向應(yīng)著重放在計算過程,以及若碰到極點的情況該如何處理方面.