基于楔形基函數(shù)的點(diǎn)插值法

童小紅, 胡 鋼

(西安理工大學(xué) 理學(xué)院, 陜西 西安 710054)

0 引言

目前，無(wú)網(wǎng)格算法求解偏微分方程數(shù)值解是近年來(lái)迅速興起的一種新型的、高效的數(shù)值方法[1-7].無(wú)網(wǎng)格算法中，相對(duì)于徑向基函數(shù)研究，楔形基函數(shù)無(wú)論在理論研究方面還是實(shí)際應(yīng)用方面都非常少[8-10].

文獻(xiàn)[11]討論了二維空間中楔形基函數(shù)插值問(wèn)題的可解性,構(gòu)造允許向量并且利用Kriging泛函的性質(zhì)給出插值問(wèn)題的誤差估計(jì).文獻(xiàn)[12]構(gòu)造了基于楔形基函數(shù)的加權(quán)最小二乘無(wú)網(wǎng)格方法，并應(yīng)用于彈性靜力學(xué)問(wèn)題.最近，文獻(xiàn)[13]討論楔形基的配點(diǎn)無(wú)網(wǎng)格法，并給出了數(shù)值解的存在惟一性.文獻(xiàn)[14,15]分別將楔形基函數(shù)無(wú)網(wǎng)格法應(yīng)用于求解對(duì)流擴(kuò)散方程.但在現(xiàn)有的文獻(xiàn)中，基于楔形基函數(shù)的點(diǎn)插值方法還沒(méi)有相關(guān)研究.考慮到點(diǎn)插值法在計(jì)算方面具有簡(jiǎn)單實(shí)用的特點(diǎn)，并針對(duì)熱傳導(dǎo)問(wèn)題，本文提出了一種新的無(wú)網(wǎng)格方法—楔形基點(diǎn)插值法.通過(guò)數(shù)值算例驗(yàn)證該算法切實(shí)可行，而且能夠達(dá)到滿意的收斂效果.

1 楔形基點(diǎn)插值無(wú)網(wǎng)格法

考慮熱傳導(dǎo)方程的邊值問(wèn)題：

(1)

其中x= (x, y)T，Ω是有界矩形區(qū)域，光滑邊界?Ω，2為L(zhǎng)aplace算子，k為常數(shù).

1.1 楔形基函數(shù)近似

定理1[1]如果φ(·)不是多項(xiàng)式，則楔形基{φ(cTx+d)}可以逼近幾乎所有的函數(shù).

(2)

(3)

根據(jù)楔形基函數(shù)的性質(zhì)，有φi(xj)=φj(xi)，所以Φ0為對(duì)稱矩陣.

1.2 數(shù)值算法的構(gòu)造

采用點(diǎn)插值法，于是由(1)式得到在時(shí)間層tn的離散：

(4)

為保證近似函數(shù)惟一性，對(duì)系數(shù)附加如下條件：

(5)

寫成矩陣形式

(6)

G為(6)式中時(shí)間層t=tn-1時(shí)右端項(xiàng)的值組成的向量.

1.3 數(shù)值解的存在惟一性

定理2[11]若楔形基函數(shù)Φ(x)的傅立葉變換F[Φ](ω)幾乎處處大于零，則插值矩陣Q是正定的.

定理3[13]函數(shù)Φ(x)=φ(a·x)是正定函數(shù)的充分必要條件是其傅里葉變換F[Φ](ω)幾乎處處大于零.

定理4若楔形基函數(shù)Φ(x)是正定的，則楔形基函數(shù)點(diǎn)插值法求解(1)式有惟一的數(shù)值解.

證明：首先證明Φ是正定的.?ζ≠0∈RN,

P的前η行構(gòu)成η級(jí)Vandermonde行列式，由于節(jié)點(diǎn)的互異性，可知前η行為其最大線性無(wú)關(guān)組，故rank(P)=η.于是可對(duì)P進(jìn)行行初等變換：P1=Q1P(Q1是一個(gè)n×n可逆矩陣)，使P1的線性無(wú)關(guān)組位于前η行，并且滿足后n-η行的元素值皆為0，有

PTΦ-1P=(Q1-1)P1TΦ-1Q1-1P1=

P1T(Q1-1)TΦ-1Q1-1

P1=P1T((Q1-1)TΦ-1Q1-1)P1

(7)

Q1可逆，故(Q1-1)TΦ-1Q1-1可逆.對(duì)P1進(jìn)行分塊：P1T=[P20]，顯然P2為η×η可逆矩陣，所以

Φ對(duì)稱正定，其合同矩陣(Q1-1)TΦ-1Q1-1對(duì)稱正定，則Φ1對(duì)稱正定，P2TΦ1-1P2對(duì)稱正定，故PTΦ-1P對(duì)稱正定.所以在時(shí)間層t=tn(5)存在惟一數(shù)值解.因此由定理3知，若楔形基函數(shù)正定， (4)存在惟一的數(shù)值解.

2 數(shù)值算例和分析

算例1 考慮一維初邊值問(wèn)題：

(8)

其中Ω={x|0≤x≤2}，解析解為u(x,t)=exp(0.001t)x(1-x).u0、g1和gN由解析解確定.采取L2誤差估計(jì).取IMQ為楔形基函數(shù).用L2范數(shù)進(jìn)行誤差估計(jì).計(jì)算方向在圓盤上取得，方向?yàn)椋篴s=(cos(θs),sin(θs)),θs=(s-1)π/9，s=1,2,…,m；時(shí)間步長(zhǎng)Δt=0.01，取m=4.計(jì)算結(jié)果如表1.

表1 本文方法數(shù)值解的誤差計(jì)算

通過(guò)表1，我們可以看到楔形基點(diǎn)插值法求解一維熱傳導(dǎo)方程能得到滿意的收斂效果.

圖1～圖4分別是本文方法、徑向基點(diǎn)插值法、有限差分法和有限元法進(jìn)行計(jì)算結(jié)果.通過(guò)圖1和圖2可以看到楔形基點(diǎn)插值法和徑向基點(diǎn)插值法一樣，都能達(dá)到滿意的收斂效果.又將計(jì)算時(shí)間和計(jì)算誤差進(jìn)行了比較，結(jié)果如表2.

表2 本文方法計(jì)算誤差與其他數(shù)值

由表2我們得到結(jié)論：本文楔形基點(diǎn)插值無(wú)網(wǎng)格法和徑向基無(wú)網(wǎng)格法都能達(dá)到滿意的收斂效果；與有限差分法和有限元法相比，本文數(shù)值算法在計(jì)算誤差方面有明顯的優(yōu)勢(shì).通過(guò)表1和表2發(fā)現(xiàn)，影響本文楔形基無(wú)網(wǎng)格點(diǎn)插值法數(shù)值解精度的因素有插入節(jié)點(diǎn)數(shù)的取法和形參數(shù)的選取.

圖1 本文的楔形基點(diǎn)插值法數(shù)值解

圖2 徑向基點(diǎn)插值法數(shù)值解

圖3 有限差分法數(shù)值解

算例2 考慮二維初邊值問(wèn)題：

(9)

其中Ω={(x,y)|0≤x≤1,0≤y≤1} ,解析解為u(x,y,t)=tsin(πx)sin(πy).u0和g由解析解確定；我們采取L2誤差估計(jì).取MQ為楔形基函數(shù).

圖4 有限元法數(shù)值解

計(jì)算方向在圓盤上取得，方向?yàn)椋篴s=(cos(θs),sin(θs)),θs=(s-1)π/(m-1), s=1,2,…,m；時(shí)間步長(zhǎng)Δt=0.01，取m=9,10,11,12,13,14,15,16.誤差計(jì)算結(jié)果如圖5～圖10、表3和表4.

圖5 本文數(shù)值解誤差t=0.5,m=9

圖6 本文數(shù)值解誤差t=1.2,m=10

圖7 本文數(shù)值解誤差t=1.8,m=12

圖8 本文數(shù)值解誤差t=3.2,m=14

圖9 本文數(shù)值解誤差t=5.5,m=14

圖10 本文數(shù)值解誤差t=16,m=16

圖11 本文數(shù)值解誤差t=0.4,m=9

tcm=9計(jì)算誤差 m=11計(jì)算誤差m=13計(jì)算誤差m=15計(jì)算誤差0.51.61.138 7e-0046.505 5e-0044.694 4e-0042.611 3e-0040.81.569.697 0e-0057.035 8e-0044.425 2e-0042.287 8e-00411.552.467 7e-0049.399 7e-0041.829 0e-0041.036 8e-0031.51.552.997 8e-0048.783 7e-0041.871 1e-0041.056 0e-00331.61.076 6e-0041.867 7e-0044.372 9e-0048.997 8e-00451.552.606 2e-0049.166 4e-0041.830 5e-0041.021 8e-003151. 54.310 5e-0042.591 4e-0041.854 0e-0041.118 5e-003

表4 本文方法數(shù)值解的誤差計(jì)算結(jié)果

由圖5～圖10我們看到，本文無(wú)網(wǎng)格方法具有很好的收斂效果.再由表3和表4的數(shù)值分析我們發(fā)現(xiàn)，影響本文楔形基無(wú)網(wǎng)格配點(diǎn)法數(shù)值解精度的因素包括方向的確定.

當(dāng)時(shí)間步長(zhǎng)Δt= 0.02，取m= 9,10,11,12.取Gaussians為楔形基函數(shù)，計(jì)算結(jié)果如表5和圖11～圖14.

圖12 本文數(shù)值解誤差t=0.4,m=10

圖13 本文數(shù)值解誤差t=0.4,m=11

圖14 本文數(shù)值解誤差t=0.4,m=12

tcm=9計(jì)算誤差 m=10計(jì)算誤差m=11計(jì)算誤差m=12計(jì)算誤差0.42.632.381 1e-0041.814 7e-0041.281 0e-0041.670 9e-0040.92.54.131 3e-0041.039 5e-0047.101 7e-0057.929 8e-0051.22.44.492 2e-0048.900 8e-0053.392 0e-0045.684 5e-0051.62.464.280 3e-0041.187 5e-0047.410 9e-0057.286 8e-00542.63.236 3e-0041.121 2e-0041.022 7e-0041.830 9e-00462.353.167 4e-0041.823 2e-0049.109 2e-0058.015 4e-005162.454.001 3e-0041.069 2e-0045.676 5e-0057.595 0e-005

由表5的數(shù)值和圖11～圖14分析我們發(fā)現(xiàn)，影響本文無(wú)網(wǎng)格法數(shù)值解精度的因素有時(shí)間步長(zhǎng)的選取和基函數(shù)的選取.

3 結(jié)論

本文所構(gòu)造的楔形基點(diǎn)插值無(wú)網(wǎng)格法求解對(duì)熱傳導(dǎo)方程切實(shí)可行，該方法是真正的無(wú)網(wǎng)格方法，并且該無(wú)網(wǎng)格法仍具有非常高的計(jì)算精度.本文無(wú)網(wǎng)格法數(shù)值解的計(jì)算精度與基函數(shù)的形參數(shù)的選取、插入節(jié)點(diǎn)的取法、時(shí)間步長(zhǎng)的選取、基函數(shù)的選取、方向的選取都相關(guān).

[1] 吳宗敏.徑向基函數(shù)、散亂數(shù)據(jù)擬合與無(wú)網(wǎng)格偏微分方程數(shù)值解[J].工程數(shù)學(xué)學(xué)報(bào),2002,19(2):1-12.

[2] 張 雄,宋康祖,陸明萬(wàn).無(wú)網(wǎng)格法研究進(jìn)展及其應(yīng)用[J].計(jì)算力學(xué)學(xué)報(bào),2003,20(6):730-742.

[3] Zhang X,Song K Z,Lu M W,et al.Meshless methods based on collocation with radial basis functions[J].Computational Mechanics,2000,26(4):333-343.

[4] Duan Y,Tan Y.Meshless collocation method based on Dirichlet-Neumann substructure iteration[J].Applied Mathematics and Computation,2005,166(2):373-384.

[5] Wang J G,Liu G R.A point interpolation meshless method based on radial basis functions[J].International Journal for Numerical Methods in Engineering,2002,54(11):1 623-1 648.

[6] 張宏偉,李美香,李衛(wèi)國(guó).關(guān)于配點(diǎn)型無(wú)網(wǎng)格法邊界條件處理技術(shù)[J].大連理工大學(xué)學(xué)報(bào),2010,50(4):614-618.

[7] 張 林,歐陽(yáng)潔,張文彬,等.基于無(wú)網(wǎng)格穩(wěn)定化方案求解非穩(wěn)態(tài)強(qiáng)對(duì)流問(wèn)題的自適應(yīng)節(jié)點(diǎn)加密技術(shù)[J].計(jì)算力學(xué)學(xué)報(bào),2010,27(3):446-450.

[8] 張立偉.有限平面波的線性組合的楔形基函數(shù)逼近[J].山東大學(xué)學(xué)報(bào)(理學(xué)版),2006,41(6):89-92.

[9] 王志剛, 秦新強(qiáng), 黨發(fā)寧, 等.楔形基無(wú)網(wǎng)格法解的存在唯一性[J].山東大學(xué)學(xué)報(bào)(理學(xué)版),2010,45(2):44-49.

[10] 田雙亮,房保言,王志剛.基于楔形基函數(shù)的微分博弈兩點(diǎn)邊值問(wèn)題的逼近[J].山東大學(xué)學(xué)報(bào)(理學(xué)版),2011,46(8):38-41.

[11] 張立偉.楔形基函數(shù)插值及其誤差估計(jì)[J].復(fù)旦學(xué)報(bào)(自然科學(xué)版),2005,44(2):302-306.

[12] 署恒木,黃朝琴,李翠偉.基于楔形基函數(shù)的一種新型無(wú)網(wǎng)格法[J].中國(guó)石油大學(xué)學(xué)報(bào),2008,32(3):109-113.

[13] Wang Z,Qin X,Wei G,et al.Meshless method with ridge basis functions[J].Applied Mathematics and Computation,2010,217(5):1 870-1 886.

[14] 童小紅, 秦新強(qiáng).一種新的求解對(duì)流擴(kuò)散邊值問(wèn)題的無(wú)網(wǎng)格方法[J].計(jì)算力學(xué)學(xué)報(bào),2012,29(5):716-720.

[15] 王麗華,秦新強(qiáng),王全九.對(duì)流占優(yōu)擴(kuò)散方程的特征線楔形基無(wú)網(wǎng)格法[J].數(shù)學(xué)雜志,2012,32(1):99-107.