Burgers方程的正交基無網(wǎng)格Galerkin解法

陳維維, 趙鳳群, 路小平

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

0 引言

無網(wǎng)格法是近年來興起的一種數(shù)值計(jì)算方法，用有限元方法解決問題時(shí)出現(xiàn)的困難，有時(shí)用無網(wǎng)格法求解則有明顯的優(yōu)勢.因此無網(wǎng)格法逐漸成為計(jì)算力學(xué)一個(gè)十分吸引學(xué)者研究的課題.它最大的優(yōu)點(diǎn)就在于部分或徹底的消除了網(wǎng)格對節(jié)點(diǎn)的約束，采用節(jié)點(diǎn)所在的影響域內(nèi)的信息構(gòu)造數(shù)值逼近.目前已經(jīng)提出了十余種無網(wǎng)格法，如：無網(wǎng)格局部Petrov-Galerkin法[1]、再生核粒子法[2]、無網(wǎng)格伽遼金法[3]、光滑粒子流體動(dòng)力法等[4].特別地，由Belytschko等人提出的無網(wǎng)格伽遼金法(Meshless Galerkin Method)已被成功的應(yīng)用到工程力學(xué)中的許多領(lǐng)域，該方法優(yōu)點(diǎn)在于采用移動(dòng)最小二乘法[5](Moving Least Squares，簡稱MLS構(gòu)造近似函數(shù))構(gòu)造的近似函數(shù)具有較高的光滑連續(xù)性，保證了計(jì)算結(jié)果不僅具有較高的精度而且還有比較好的穩(wěn)定性，缺點(diǎn)是在計(jì)算形函數(shù)和形函數(shù)的導(dǎo)數(shù)過程中都要計(jì)算矩陣的逆，這樣又會(huì)導(dǎo)致計(jì)算效率低，當(dāng)基函數(shù)的項(xiàng)數(shù)較大時(shí)容易使方程組產(chǎn)生病態(tài).為了克服無網(wǎng)格Galerkin法的上述弊端，本文在無網(wǎng)格伽遼金法中使用正交基函數(shù)，使矩陣求逆變得簡單，不但避免了產(chǎn)生病態(tài)方程組的可能，而且提高了計(jì)算效率.對于Burgers方程，在時(shí)間域上用θ加權(quán)法[6]進(jìn)行離散，而空間域上采用正交基無網(wǎng)格Galerkin法進(jìn)行離散，構(gòu)造了一種新的數(shù)值計(jì)算方法θ加權(quán)-正交基無網(wǎng)格Galerkin法.

1 無網(wǎng)格Galerkin法的基本理論

1.1 移動(dòng)最小二乘近似

在求解域Ω中，函數(shù)u(x)可近似為uh(x)，求解域中插值節(jié)點(diǎn)的集合為{x1,x2,…,xN}，則函數(shù)u(x)的移動(dòng)最小二乘近似式uh(x)可以定義為[7,8]

(1)

式中uh(x)為函數(shù)u(x)的近似表達(dá)式.pi(x)基函數(shù)，m是基函數(shù)的個(gè)數(shù)，系數(shù)ai(x),i=1,2,…,m是空間坐標(biāo)x的函數(shù).

為了求得系數(shù)a(x)，我們可以采用局部近似的加權(quán)最小二乘法擬合求解，即使局部近似誤差泛函J取極小值

(2)

式中，w(x-xi)是帶有緊支撐性質(zhì)的光滑連續(xù)權(quán)函數(shù)，即在x的影響域內(nèi)部，w(x-xi)>0，在其邊界和外部，w(x-xi)=0，n是積分點(diǎn)x鄰域內(nèi)的節(jié)點(diǎn)數(shù).

由于J取極小值，由極值原理可得

(3)

其中的矩陣A(x)、B(x)分別是：

A(x)=pTw(x)p′B(x)=pTw(x)

當(dāng)矩陣A為非奇異矩陣時(shí)，由(3)可解得

a(x)=A-1(x)B(x)u

(4)

將a(x)代入(1)式中，得到逼近函數(shù)uh(x)的表達(dá)式是

uh(x)=pT(x)A-1(x)B(x)u=φ(x)u

(5)

其中，φ(x)=pT(x)A-1(x)B(x)，φ(x)是MLS形函數(shù).

1.2 基于正交基函數(shù)的移動(dòng)最小二乘近似

在上述移動(dòng)最小二乘法中，如果基函數(shù)的個(gè)數(shù)m較大時(shí)，方程(3)有可能產(chǎn)生病態(tài)，甚至有可能是奇異的，這樣就導(dǎo)致該方程很難求出解.為了克服這一問題，可將基函數(shù)取為正交函數(shù)，這樣不僅可以避免所得方程病態(tài)或者奇異，而且避免了求矩陣的逆，簡化了計(jì)算.

對于權(quán)函數(shù){wi}和點(diǎn)集{xi}，如果函數(shù)p1(x),p2(x),…,pm(x)滿足條件

(6)

如任意給出一組基函pk(x) (k=1,2,…,m) 通過下述過程可將基函數(shù)pk(x)轉(zhuǎn)化為正交基函數(shù)qk(x)

令q1(x)=p1(x)

?

(7)

由此(4)式中的矩陣A可化為

(8)

這時(shí)將對角化后的矩陣A帶入(4)式可求解出系數(shù)

(9)

將系數(shù)aj(x)代入(5)式可得

uh(x) =qT(x)A-1(x)B(x)ui

(10)

其中Ni(x)是構(gòu)造好的新的形函數(shù).

為了驗(yàn)證正交基函數(shù)的移動(dòng)最小二乘近似的有效性，取如下的函數(shù)進(jìn)行計(jì)算分析：

f(x,y)=(2x+sin(x))3·(ycos(y))4

本文選取高斯函數(shù)作為權(quán)函數(shù)，基函數(shù)取pT(x,y)=(1,x,y)，使用20×20個(gè)均勻分布節(jié)點(diǎn)進(jìn)行計(jì)算，圖1和圖2分別給出了上述函數(shù)關(guān)于x和y的偏導(dǎo)的近似結(jié)果.

圖1 f關(guān)于x的偏導(dǎo)數(shù)逼近(y=11)

圖2 f關(guān)于y的偏導(dǎo)數(shù)逼近(x=8)

由圖1和圖2可以看出，采用正交基函數(shù)的移動(dòng)最小二乘近似方法進(jìn)行函數(shù)逼近時(shí)，函數(shù)偏導(dǎo)數(shù)的數(shù)值解和精確解逼近效果很好，這就證明了該方法能很好地應(yīng)用于函數(shù)逼近的問題中，同時(shí)也表明了該方法的有效性.

2 偏微分方程的無網(wǎng)格Galerkin法離散

考慮如下形式的Burgers方程

(11)

首先，采用θ加權(quán)法對方程(11)的時(shí)間域進(jìn)行離散，則有

(12)

將(12)式代入方程(11)則得到關(guān)于時(shí)間域離散的半離散方程

(13)

其次，采用θ加權(quán)-正交基無網(wǎng)格Galerkin法對空間域進(jìn)行離散，設(shè)uh(x)=Nui，并用罰函數(shù)處理邊界條件，則原方程可化為如下方程組：

(k1-Δtθk2)un+1=[Δt(1-θ)k2+Δtk3-k1]un

(14)

3 數(shù)值算例

圖3 不同時(shí)刻的數(shù)值解

圖3給出了當(dāng)t=0.0,0.01,0.1,0.15,0.2,0.25時(shí)用本文方法計(jì)算出的結(jié)果.由圖3可以看出，隨著時(shí)間步長的變化，用θ加權(quán)-正交基無網(wǎng)格Galerkin法求解Burgers方程，可以有效的避免數(shù)值振蕩.

表1給出了粘性系數(shù)ε=0.01在不同點(diǎn)處不同時(shí)刻的數(shù)值解，并與精確解和文獻(xiàn)[9]的算法數(shù)值解進(jìn)行了比較，可以看出本文的數(shù)值解和精確解很吻合，表明本文的方法是很有效的.

表1 ε=0.01,N=40,Δt=0.000 2時(shí)刻所得的數(shù)值解及與文獻(xiàn)[9]的對比

4 結(jié)束語

本文在傳統(tǒng)的無網(wǎng)格Galerkin法中引入了正交基函數(shù)，形成正交基函數(shù)的無網(wǎng)格伽遼金法，對時(shí)間域的離散引入了θ加權(quán)法，構(gòu)造出了一種θ加權(quán)-正交基函數(shù)無網(wǎng)格伽遼金法，并成功的將該方法應(yīng)用于求解Burgers方程中.由于Burgers方程是典型的拋物型偏微分方程，所以θ加權(quán)-正交基函數(shù)的無網(wǎng)格Galerkin法也可以作為一種有效的數(shù)值方法用于求解其他的拋物型偏微分方程.此外，由于正交基函數(shù)的無網(wǎng)格方法簡化了數(shù)值計(jì)算的前處理過程，所以選用正交基函數(shù)的無網(wǎng)格Galerkin方法求解一些偏微分方程也具有一定的優(yōu)勢.

[1] S.N.Atluri,T.Zhu.A new meshless local petrov-galerkin(MLPG) approach in computa-tion mechanics[J].Computational Mechanics,1998,22(2):117-127.

[2] W.K.Liu,S.Jun,Y.F.Zhang.Reproducing kernel particale methods[J].Numerical Method in Fluids,1995,20(8):1 081-1 106.

[3] T.Belytschko,Y.Y.Lu,L.Gu.Element-free Galerkin method[J].Numerical Methods in Engineering,1994,37(2):229-256.

[4] M.R.Bate,A.Burkert.Resolution requiements for smoothed particle hydrodynamics calculations with self-gravity[J]. Monthly Notices of the Royal Astronomical Society,1997,228(4):1 060-1 072.

[5] 張 雄，劉 巖.無網(wǎng)格法[M].北京：清華大學(xué)出版社，2004.

[6] J.Done,A.Hnerta.Finite Element Methods for Flow Problems[M].England:Wiley,2003.

[7] Belytschko T,Krongauz Y,Organ D.Meshless methods:an over view and recent developments[J].Computer Methods in Applied Mechanics and Engineering,1996,139(1):3-47.

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

[9] 劉萬海，孫健安.用五次B樣條Galerkin有限元方法求Burgers方程的數(shù)值解[J].西北師范大學(xué)學(xué)報(bào)，2009，45(2)：35-38.