朱雙彪

(南京財(cái)經(jīng)大學(xué)應(yīng)用數(shù)學(xué)學(xué)院,江蘇 南京 210023)

許多科學(xué)問(wèn)題都可以轉(zhuǎn)化為偏微分方程問(wèn)題.通過(guò)數(shù)值求解這類(lèi)問(wèn)題的方法有很多,自然邊界元[1-3]與有限元耦合法就是其中的一類(lèi).耦合法的好處在于,2種方法結(jié)合時(shí)能夠彌補(bǔ)各自的不足.利用耦合法不僅可以求解圓和橢圓區(qū)域問(wèn)題,還可以求解凹角區(qū)域問(wèn)題[4-5].但在求解耦合問(wèn)題時(shí),由于是用直邊三角形作有限元分析的,因此會(huì)產(chǎn)生誤差.采用曲邊有限元即曲邊三角形代替直邊三角形[6-8],就可以大大降低誤差.筆者擬給出凹角區(qū)域上求解泊松方程邊值問(wèn)題的曲邊有限元和自然邊界元的耦合法.

1 曲邊有限元與自然邊界元耦合法

如圖1所示,取夾角α、線段Γ0和Γ1、光滑的圓弧Γ，一起圍成區(qū)域Ω.在區(qū)域Ω內(nèi)求解如下邊值問(wèn)題：

圖1 區(qū)域ΩFig. 1 Area Ω

(1)

(2)

用半徑為R的人工圓弧Γ2置于區(qū)域Ω內(nèi),其中Γ2={(R,θ)|0≤θ≤α},且dist(Γ,Γ2)>0.Γ2將區(qū)域Ω分為2個(gè)區(qū)域,Γ2包圍的區(qū)域?yàn)棣?，Γ2和Γ包圍的區(qū)域?yàn)棣?,如圖2所示.在區(qū)域Ω2內(nèi)用自然邊界元方法,在區(qū)域Ω1內(nèi)用曲邊有限元方法.

圖2 區(qū)域Ω引入人工邊界Fig. 2 Artificial Boundary of Area Ω

在耦合法中，一般是用直邊三角形近似劃分有限區(qū)域,但這種處理方式存在一定的誤差.采用曲邊元素,就可以很好地?cái)M合待求解區(qū)域彎曲的邊界[4-5].

(3)

取參考單元G,其頂點(diǎn)g1=(0,0),g2=(1,0),g3=(0,1).記面積坐標(biāo)(λ1,λ2,λ3),于是得到如下線性變換：

其中λ1=1-ξ-η,λ2=ξ,λ3=η.作如下變換F：

(4)

其中

(5)

(6)

這里

線性變換(3)～(6)可以將參考單元的3條邊一一映射為曲邊三角形的2條直邊和1個(gè)曲邊[6].

σ={Ni:i=1,2,3}，是一組線性函數(shù).定義Ni(φ)=φ(pi),φ∈Σ,i=1,2,3.在τ和有限元正則的情況下,各個(gè)曲邊三角單元的插值誤差估計(jì)為

(7)

2 曲邊有限元與自然邊界元凹角耦合法的收斂性和誤差估計(jì)

變分問(wèn)題(2)的離散問(wèn)題是

(8)

變分問(wèn)題(8)有唯一的連續(xù)依賴(lài)于f的解.

引理1[1]

證明因?yàn)镾h∈H1(Ω1),所以在

中取v=vh∈Sh,再與

相減,得到

又因?yàn)閡h∈Sh,所以

證畢.

由引理1和(7)式可以得到定理1,2：

定理1(收斂性) 假設(shè)插值算子滿(mǎn)足

那么

定理2(能量模估計(jì)) 設(shè)u∈H2(Ω1),那么

由定理1給出的近似解滿(mǎn)足收斂性可知，通過(guò)自然邊界元和曲邊有限元耦合法求解問(wèn)題(1)是可行性.

‖u-uh‖L2(Ω1)≤Ch2‖u‖2,Ω1.

對(duì)比能量模估計(jì)和L2模估計(jì)，L2模估計(jì)是最優(yōu)的.

3 數(shù)值實(shí)例

例1在圖1中取夾角α=7π/4.采用曲邊有限元與自然邊界元耦合法求解如下外邊值問(wèn)題：

其中Γ0,Γ1,Γ的邊界條件分別為

Γ0={(r,θ)|r≥3,θ=0},Γ1={(r,θ)|r≥3,θ=α},Γ={(r,θ)|r=3,0<θ<α},

人工邊界為

Γ2={(r,θ)|r=1.5,0<θ<α}.

通過(guò)Matlab軟件對(duì)區(qū)域進(jìn)行加密剖分，結(jié)果如圖3，4所示.圖3是將邊界劃分為33個(gè)節(jié)點(diǎn)，圖4是將邊界劃分為17個(gè)節(jié)點(diǎn)，可以看出圖3的剖分比圖4更密.

圖3 mesh=8×33Fig. 3 mesh=8×33

圖4 mesh=4×17Fig. 4 mesh=4×17

取剖分網(wǎng)格N=17,33,65,129.有限元采用直邊三角元和曲邊三角元產(chǎn)生的誤差分別見(jiàn)表1,2.

表1 有限元采用直邊三角元產(chǎn)生的誤差Table 1 Error Value of Straight Triangular Element in Finite Element

表2 有限元采用曲邊三角元產(chǎn)生的誤差Table 2 Error Value of Curved Edge Element in Finite Element

從表1,2可知：誤差隨著剖分網(wǎng)格的加密越來(lái)越??；當(dāng)剖分網(wǎng)格相同時(shí),L2誤差階都接近于ο(h2),但有限元采用曲邊三角元產(chǎn)生的誤差比直邊三角元的小.這說(shuō)明，曲邊有限元采用曲邊三角單元代替原來(lái)的直邊三角形，可以克服邊界近似產(chǎn)生的誤差，降低數(shù)值積分誤差.