夏茂輝，趙玉鳳，呂 鵬，翟育鵬，任偉和

(燕山大學(xué)理學(xué)院，河北秦皇島066004)

0 引言

近年來(lái)，無(wú)網(wǎng)格法已經(jīng)發(fā)展成為一種功能較強(qiáng)、穩(wěn)定性高的數(shù)值計(jì)算方法，而無(wú)網(wǎng)格迦遼金法(EFGM)是發(fā)展最為成熟的一種無(wú)網(wǎng)格法.傳統(tǒng)的移動(dòng)最小二乘法構(gòu)造出一個(gè)近似函數(shù)［1］，它的光滑性好且導(dǎo)數(shù)連續(xù)，其形函數(shù)及其導(dǎo)數(shù)僅與節(jié)點(diǎn)有關(guān)，但此方法容易導(dǎo)致系統(tǒng)方程病態(tài)，于是在此近似方案中采用帶權(quán)的正交基函數(shù)，形成一種改進(jìn)的移動(dòng)最小二乘近似(IMLS).IMLS近似與EFGM相結(jié)合構(gòu)成了IEFGM，該方法避免了矩陣求逆，且計(jì)算量小、數(shù)值穩(wěn)定性好，不會(huì)形成病態(tài)方程組，因此有更高的計(jì)算速度和精度，在彈性力學(xué)和斷裂力學(xué)中已得到了充分應(yīng)用［2-5］，但在溫度場(chǎng)中還未得到深入研究.筆者在文獻(xiàn)［6］的基礎(chǔ)上對(duì)移動(dòng)最小二乘近似進(jìn)行改進(jìn)，用IEFG法求解穩(wěn)態(tài)熱傳導(dǎo)問(wèn)題，并將其方法應(yīng)用到矩形域和圓形域中，并與EFG法和精確解［7-8］進(jìn)行比較，應(yīng)用MATLAB編程將其結(jié)果進(jìn)行可視化處理，表明該方法是一種收斂快、精度高、簡(jiǎn)便有效的方法.

1 移動(dòng)最小二乘近似［1］

考慮點(diǎn)x的鄰域Ωx，它位于全域Ω內(nèi)，為求近似函數(shù)，在子域Ωx中隨機(jī)分布有限個(gè)節(jié)點(diǎn)xi(i=1，2，…，n)，函數(shù) T 的近似表達(dá)式Th(x)，?x∈Ωx，可以定義為:

式中:PT(x)=［p1(x)p2(x)… pm(x)］是m次完備單項(xiàng)式的基;a(x)是包含系數(shù)ai(x)，i=1，2，…，m的向量，這些系數(shù)是空間坐標(biāo)x=［x y］的函數(shù).對(duì)于二維問(wèn)題:

線性基

二次基

定義加權(quán)離散L2模為:

式中:wi(x)是節(jié)點(diǎn)i在x的權(quán)函數(shù);Ti(i=1，2，…，n)是節(jié)點(diǎn)i處的名義節(jié)點(diǎn)值，一般它不是未知試函數(shù)Th(x)的節(jié)點(diǎn)值;n為插值基點(diǎn)數(shù).

式(2)可用矩陣形式表達(dá)為:

其中，

2 改進(jìn)的移動(dòng)最小二乘近似

在MLS近似中，式(4)有時(shí)是病態(tài)的，甚至出現(xiàn)奇異現(xiàn)象，而且難以得到數(shù)值解.為克服此缺點(diǎn)，采用 IMLS［2-5］近似.

對(duì)于?f(x)，g(x)∈span(p)，定義

其中，span(p)是Hilbert空間;(f，g)是內(nèi)積.對(duì)于一系列點(diǎn){xi}和權(quán)函數(shù){wi}，若函數(shù)p1(x)，p2(x)…pn(x)∈span(p)且滿(mǎn)足:

所以由式(4)、式(7)可得:

其中，I=1，2，…，n.若 pi(x)∈span(p)，i=1，2，…，m;(pi，pj)=0，i≠j，則(9)式為:

于是可得ai( x)為:

3 改進(jìn)型EFG法(IEFGM)

考慮二維穩(wěn)態(tài)熱傳導(dǎo)問(wèn)題.

式中:T表示溫度;T0為初始溫度;Γφ為狄利克雷(Dirchlet)邊界;Γq為諾依曼(Neumann)邊界;ρ為材料密度;kx、ky是材料沿x、y方向的熱傳導(dǎo)系數(shù);Q為物體內(nèi)部的熱源密度;nx和ny是邊界外法線的方向余弦;ˉq是給定的熱流量.

與方程(17)相對(duì)應(yīng)的弱形式為

把式(14)代入式(18)，可得到離散方程:

由IMLS可得節(jié)點(diǎn)溫度值

4 數(shù)值算例

算例1 圓域上的穩(wěn)態(tài)熱傳導(dǎo)問(wèn)題.考慮單位圓域上的穩(wěn)態(tài)溫度分布場(chǎng)，其控制方程為:

邊界條件為:

精確解的直角坐標(biāo)表達(dá)式為［7］:

如圖1在圓域內(nèi)布置了120個(gè)規(guī)則節(jié)點(diǎn)，高斯積分背景網(wǎng)格與節(jié)點(diǎn)相重合，每個(gè)網(wǎng)格采用2×2高斯積分，采用三次樣條權(quán)函數(shù)和帶權(quán)的線性正交基函數(shù).由IEFG法得直角坐標(biāo)系下φ=0半徑線段上溫度曲線對(duì)比如圖2，溫度分布如圖3，將改進(jìn)型無(wú)網(wǎng)格迦遼金法(IEFGM)和精確解進(jìn)行比較.

算例2 三維矩形域上的穩(wěn)態(tài)熱傳導(dǎo)問(wèn)題.立方體邊長(zhǎng)a=b=c=1，其控制方程和邊界條件如下［8］:

采用均勻布點(diǎn)，布置945個(gè)節(jié)點(diǎn)對(duì)立方體進(jìn)行離散(165個(gè)背景網(wǎng)格)，影響域取為立方體影響域，每個(gè)背景網(wǎng)格內(nèi)采用64個(gè)高斯積分點(diǎn)，權(quán)函數(shù)為三次樣條權(quán)函數(shù)，二次正交基函數(shù).圖4為z=1時(shí)的溫度分布對(duì)比圖.圖5為z=1時(shí)等溫線對(duì)比圖.將改進(jìn)型無(wú)網(wǎng)格迦遼金法(IEFGM)和精確解進(jìn)行比較.

由以上數(shù)值算例可以看出，通過(guò)IEFG法求解穩(wěn)態(tài)熱傳導(dǎo)問(wèn)題，結(jié)果比無(wú)網(wǎng)格迦遼金法精度高，和精確解高度吻合.得知IEFG法可以有效地用于求解矩形域和圓域內(nèi)穩(wěn)態(tài)熱傳導(dǎo)問(wèn)題.

5 結(jié)論

(1)無(wú)網(wǎng)格迦遼金法僅需要布置節(jié)點(diǎn)，不需要?jiǎng)澐謫卧?，其結(jié)果與有限元高度吻合，而改進(jìn)型無(wú)網(wǎng)格迦遼金法比迦遼金法有更高的精度和效率，且不會(huì)導(dǎo)致系統(tǒng)方程產(chǎn)生病態(tài).

(2)嘗試將改進(jìn)的無(wú)網(wǎng)格迦遼金法應(yīng)用到求解矩形域和圓域穩(wěn)態(tài)熱傳導(dǎo)問(wèn)題，有較好的精度和效率，且結(jié)果與精確解［7-8］高度吻合，在溫度場(chǎng)中具有廣泛的應(yīng)用前景.

［1］ 周暉，李勇.無(wú)網(wǎng)格伽遼金法在二維土體沉降分析中的應(yīng)用［J］.湖南科技大學(xué)學(xué)報(bào):自然科學(xué)版，2009，24(2):53-56.

［2］ 程玉民.移動(dòng)最小二乘法研究進(jìn)展與評(píng)述［J］.計(jì)算機(jī)輔助工程，2009，18(2):5-11.

［3］ ZHANG Zan，LIEW K M，CHENG Yu-min.Coupling of the improved element-free Galerkin and boundary element methods for two-dimensional elasticity problems［J］.Eng Anal Boundary Elem，2008(32):100-107.

［4］ 邊燕飛.改進(jìn)型無(wú)網(wǎng)格伽遼金法(IEFG)的研究及其應(yīng)用［J］.合肥工業(yè)大學(xué)學(xué)報(bào):自然科學(xué)版，2009，32(4):32-35.

［5］ ZHANG Zan，LIEW K M，CHENG Yu-min.Analyzing 2D fracture problems with the improved elementfree Galerkin method［J］.Eng Anal Boundary Elem，2008，32:241-250.

［6］ 高會(huì)生.MATLAB原理與工程應(yīng)用［M］.北京:電子工業(yè)出版社，2005.

［7］ 商立英.非結(jié)構(gòu)化網(wǎng)格上非穩(wěn)態(tài)導(dǎo)熱問(wèn)題的求解與應(yīng)用［D］.南京:南京理工大學(xué)動(dòng)力工程學(xué)院，2006.

［8］ 陳琳.無(wú)網(wǎng)格在熱傳導(dǎo)問(wèn)題中的應(yīng)用［D］.南京:南京理工大學(xué)動(dòng)力工程學(xué)院，2009.