Helmholtz方程的楔形基區(qū)域分解法

童小紅，秦新強(qiáng)

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

Helmholtz方程的楔形基區(qū)域分解法

童小紅，秦新強(qiáng)

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

1 引言

Helmholtz型方程在許多工程領(lǐng)域都會(huì)遇到，如電磁場(chǎng)中的波導(dǎo)問(wèn)題、薄膜振動(dòng)問(wèn)題以及海洋工程中水波衍射問(wèn)題等都是由Helmholtz型方程控制的[1-2]。因此，構(gòu)造穩(wěn)定、快速實(shí)用的數(shù)值算法求解Helmholtz型方程，有著重要的理論和實(shí)際意義。

考慮二維Helmholtz型方程的邊值問(wèn)題：

其中，Ω是有界矩形區(qū)域，光滑邊界?Ω，Δ為L(zhǎng)aplace算子，k為波數(shù)。

目前，求解Helmholtz型方程的數(shù)值方法有有限元法[3]、有限差分[4]、光譜法[5]。而無(wú)網(wǎng)格法是近年來(lái)迅速興起的一種新型、高效的數(shù)值方法，基于無(wú)網(wǎng)格方法求解偏微分方程的研究也有了一些成果[6-7]?，F(xiàn)在國(guó)際上常用的處理多元問(wèn)題的函數(shù)基有兩種：徑向基函數(shù)和楔形基函數(shù)?；趶较蚧瘮?shù)插值求偏微分方程數(shù)值解的研究已經(jīng)有了一系列成果[8-9]，而對(duì)于楔形基函數(shù)的研究，無(wú)論在理論研究方面還是實(shí)際應(yīng)用方面都非常少[10-14]。

但是，使用無(wú)網(wǎng)格配點(diǎn)法求解區(qū)域較大的問(wèn)題，當(dāng)空間步長(zhǎng)較小時(shí)，節(jié)點(diǎn)數(shù)的增加會(huì)引起系數(shù)矩陣的條件數(shù)太大，以至于出現(xiàn)病態(tài)問(wèn)題[10，15]。針對(duì)這種病態(tài)問(wèn)題，本文引入了楔形基無(wú)網(wǎng)格配點(diǎn)的區(qū)域分解方法，并應(yīng)用于求解Helmholtz型方程數(shù)值解問(wèn)題。通過(guò)數(shù)值算例表明該算法兼有楔形基無(wú)網(wǎng)格配點(diǎn)法的優(yōu)點(diǎn)，降低了系數(shù)矩陣的條件數(shù)，提高了算法的穩(wěn)定性，并達(dá)到了滿(mǎn)意的收斂效果。

2 楔形基配點(diǎn)無(wú)網(wǎng)格法

2.1 數(shù)值方法的構(gòu)造

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

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

由定理1，設(shè)u(x)的近似值u～(x)可以表示為：

其中，Ii∈W，Bi∈?Ω，λi為未知系數(shù)，Ni表示域內(nèi)的節(jié)點(diǎn)數(shù)，Nb表示邊界上的節(jié)點(diǎn)數(shù)，as為固定方向[10]，m為圓盤(pán)方向的總組數(shù)[10]，?為楔形基函數(shù)。

下面，利用楔形基函數(shù)作為逼近函數(shù)構(gòu)造無(wú)網(wǎng)格方法對(duì)問(wèn)題（1）進(jìn)行數(shù)值求解，離散方案采用配點(diǎn)法。要求式（2）在邊界處滿(mǎn)足邊界條件，在域內(nèi)滿(mǎn)足方程，由（1）可以得到：

當(dāng)j=1，2，…，(Nb+Ni)時(shí)，式（3）和（4）表示的是線性方程組，也就是橢圓型方程的楔形基配點(diǎn)無(wú)網(wǎng)格計(jì)算方法。為討論該方法解的存在唯一性，定義矩陣：

其中，A是(Nb+Ni)×(Nb+Ni)矩陣，Λ是待定系數(shù)向量，B是給定的右端向量。

2.2 數(shù)值解的存在唯一性

定理3[10]函數(shù)Ф(x)=?(a·x)是正定函數(shù)的充分必要條件是其傅里葉變換F[Ф](w)幾乎處處大于零。

定理4如果楔形基函數(shù)Ф(x)=?(a·x)是正定的，并且存在，那么AΛ=B存在唯一解。

證明首先證明和存在性。

一般地，取楔形基函數(shù)Ф(x)=?(a·x)為偶函數(shù)，由定理1和定理3可知，A1是正定的，即存在。

式（6）兩邊取行列式得：

3 楔形基函數(shù)配點(diǎn)區(qū)域分解法

假設(shè)Ω被分解為兩個(gè)重疊的區(qū)域Ω1和Ω2，記Γ=Ω1∩Ω2≠?，?2=?Ω1∩Ω2，?1=?Ω2∩Ω1，如圖1所示。

圖1 重疊區(qū)域圖

引理1[16]如果對(duì)域Ω1和Ω2的Dirichlet問(wèn)題皆有解，則對(duì)Ω的問(wèn)題也有解。

由定理4和引理1，得到下面定理：

定理5如果楔形基函數(shù)Ф(x)=?(a·x)是正定的，并且存在，那么楔形基函數(shù)配點(diǎn)區(qū)域分解算法求解Helmholtz型方程（1）存在唯一解。

下面給出求解Helmholtz型方程的楔形基函數(shù)配點(diǎn)區(qū)域分解算法：

步驟1選初始∈C(?1)，n:=0。

步驟2在Ω1上采用楔形基函數(shù)配點(diǎn)區(qū)域分解求解Helmholtz型方程的Dirichlet問(wèn)題，解得，即

然后，在Ω2上采用楔形基函數(shù)配點(diǎn)區(qū)域分解求解Helmholtz型方程的Dirichlet問(wèn)題，解得u～2n+2∈H1(Ω2)，即

步驟3，預(yù)先設(shè)定的精度為ε0，當(dāng)εn≤ε0時(shí)，結(jié)束程序；否則轉(zhuǎn)步驟2。

步驟4置n:=n+1轉(zhuǎn)步驟2。

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

算例考慮二維Helmholtz型方程初邊值問(wèn)題：

為了比較楔形基無(wú)網(wǎng)格算法的優(yōu)劣，分別用徑向基配點(diǎn)法、有限差分法和有限元法進(jìn)行計(jì)算，結(jié)果如圖2～5所示。首先，通過(guò)圖2和圖3可以看到楔形基配點(diǎn)法和徑向基配點(diǎn)法一樣，都能達(dá)到滿(mǎn)意的收斂效果，其次，將計(jì)算時(shí)間和計(jì)算誤差進(jìn)行了比較，結(jié)果如表1。

表1 本文方法與其他數(shù)值方法的比較（節(jié)點(diǎn)數(shù)13×13）

通過(guò)表1可以得到結(jié)論：楔形基配點(diǎn)法和徑向基配點(diǎn)法都能達(dá)到較高的精度，得到比較滿(mǎn)意的結(jié)果；相對(duì)于傳統(tǒng)的有限差分法和有限元法，本文的楔形基配點(diǎn)無(wú)網(wǎng)格方法更有些優(yōu)勢(shì)。

將區(qū)域Ω分解為Ω1={(x，y)|0≤x≤3/4，0≤y≤1}和Ω2={(x，y)|1/4≤x≤1，0≤y≤1}兩個(gè)區(qū)域，?1為邊x=3/4，?2為邊x=1/4。當(dāng)ε0=10-6，計(jì)算結(jié)果如圖6～9和表2所示。

通過(guò)圖6～9和表2，可以看到楔形基配點(diǎn)區(qū)域分解法求解Helmholtz型方程切實(shí)可行。再通過(guò)表2比較發(fā)現(xiàn)：楔形基配點(diǎn)區(qū)域分解法降低了系數(shù)矩陣的條件數(shù)，同時(shí)降低了計(jì)算誤差，并達(dá)到了滿(mǎn)意的收斂效果。

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

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

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

圖6 本文區(qū)域分解左區(qū)域數(shù)值解

圖7 左區(qū)域的誤差圖

圖8 本文區(qū)域分解右區(qū)域數(shù)值解

圖9 右區(qū)域的誤差圖

表2 本文數(shù)值方法的計(jì)算誤差和系數(shù)矩陣條件數(shù)

5 結(jié)束語(yǔ)

楔形基配點(diǎn)區(qū)域分解法求解Helmholtz型方程切實(shí)可行，通過(guò)數(shù)值算例和分析發(fā)現(xiàn)：基于楔形基函數(shù)配點(diǎn)區(qū)域分解法降低了系數(shù)矩陣的條件數(shù)，同時(shí)降低了計(jì)算誤差，并達(dá)到滿(mǎn)意的收斂效果。數(shù)值結(jié)果表明本文算法具有可行性，能有效地應(yīng)用于大型區(qū)域的工程問(wèn)題中。

[1]Jin B，Marin L.The plane wave method for inverse problems associated with Helmholtz-type equations[J].Engineering Analysis with Boundary Elements，2008，32（3）：223-240.

[2]Cheng H，Huang J，Leiterman T.An adaptive fast solver for the modifiedHelmholtzequationintwodimensions[J].Journal of Computational Physics，2006，211（2）：616-637.

[3]Van T，Wood A.A time-domain finite element method for Helmholtzequations[J].JournalofComputationalPhysics，2002，183（2）：486-507.

[4]Lambe L A，Luczak R，Nehrbass J W.A new finite difference method for the Helmholtz equation using symbolic computation[J].International Journal of Computational Engineering Science，2003，4（1）：121-144.

[5]Auteri F，Quartapelle L.Galerkin-legendre spectral method for the 3D Helmholtz equation[J].Journal of Computational Physics，2000，161（2）：454-483.

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

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

[8]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.

[9]Duan Y，Tan Y.Meshless collocation method based on Dirichlet-Neumann substructure iteration[J].Applied Mathematics and Computation，2005，166（2）：373-384.

[10]Wang Z，Qin X，Wei G，et al.Meshless method with ridge basisfunctions[J].AppliedMathematicsandComputation，2010，217（5）：1870-1886.

[11]Tong X，Qin X.A new meshless method for convection diffusion boundary value problems[J].Chinese Journal of Computational Mechanics，2012，29（5）：716-720.

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

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

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

[15]房保言，王志剛，梁波，等.修正Helmholtz方程的小波配點(diǎn)區(qū)域分解法[J].武漢理工大學(xué)學(xué)報(bào)，2010，32（3）：120-123.

[16]呂濤，石濟(jì)民，林振寶.區(qū)域分解算法-偏微分方程數(shù)值解新技術(shù)[M].北京：科學(xué)出版社，1992.

TONG Xiaohong,QIN Xinqiang

School of Sciences,Xi’an University of Technology,Xi’an 710054,China

A meshless collocation method for Helmholtz-type equations is developed,which is based on domain decomposition method and ridge basis function.This method overcomes the ill-conditioned problems caused by much larger condition number of the coefficient matrix when collocation method is used for solving Helmholtz-type equations defined large-scale regions. According to numerical example and analysis,it is seen that the method Helmholtz-type equations slows down the condition number of the coefficient matrix and calculation errors,and can achieve satisfactory results.

Helmholtz-type equations;collocation method;ridge basis function;domain decomposition method

基于楔形基函數(shù)和無(wú)網(wǎng)格配點(diǎn)法，提出了一種求解Helmholtz型方程區(qū)域分解法。該方法克服了在求解大規(guī)模問(wèn)題時(shí)用一般的全域配點(diǎn)法所帶來(lái)的配置矩陣為非對(duì)稱(chēng)滿(mǎn)陣，且高度病態(tài)的問(wèn)題。通過(guò)數(shù)值結(jié)果表明，該算法在求解Helmholtz型方程降低系數(shù)矩陣條件數(shù)的同時(shí)，也能夠降低誤差，并達(dá)到滿(mǎn)意的收斂效果。

Helmholtz型方程；配點(diǎn)法；楔形基函數(shù)；區(qū)域分解法

A

O225

10.3778/j.issn.1002-8331.1210-0240

TONG Xiaohong,QIN Xinqiang.Domain decomposition method for Helmholtz equations based on ridge basis function. Computer Engineering and Applications,2013,49（13）：40-42.

陜西省自然科學(xué)基金（No.2012JM1008）。

童小紅（1971—），男，講師，研究方向：計(jì)算機(jī)圖形學(xué)，偏微分方程數(shù)值解。E-mail：xiaohtong@163.com

2012-10-24

2012-12-20

1002-8331（2013）13-0040-03

CNKI出版日期：2012-12-26http://www.cnki.net/kcms/detail/11.2127.TP.20121226.1120.006.html