許文政,張耀明
(山東理工大學(xué) 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院,山東 淄博 255049)
邊界元法是分析科學(xué)和工程應(yīng)用問題的強(qiáng)有力的數(shù)值模擬工具,廣泛應(yīng)用于位勢流,彈性靜、動(dòng)力學(xué),聲傳播等線性問題及各種非均勻材料與非線性問題[1-6].但是,不同于域型方法,邊界元法需要計(jì)算奇異核積分,它的有效計(jì)算是邊界法成功實(shí)施的關(guān)鍵.時(shí)至今日,這方面的研究仍在進(jìn)行中[1-3].
現(xiàn)有許多有效計(jì)算奇異積分的方法[7-9].大致可分為直接計(jì)算方法和間接計(jì)算方法兩大類[10-20].其中,Guiggiani等[10]提出的方法具有代表性,該方法在內(nèi)蘊(yùn)坐標(biāo)系的參數(shù)平面內(nèi)操作,能計(jì)算各種階奇異積分.該方法需要將被積函數(shù)中的每個(gè)量展成單元局部距離Ω的泰勒級數(shù),因此需要許多復(fù)雜的推導(dǎo)工作和與單元有關(guān)的計(jì)算.高效偉[11]也提出了一種計(jì)算奇異積分的一般性算法,它將積分核中的非奇異部分表示成在內(nèi)蘊(yùn)坐標(biāo)系統(tǒng)下的局部距離參數(shù)的冪級數(shù),然后解析計(jì)算.確定冪級數(shù)的系數(shù)時(shí),需要求解線性方程組,而且這樣的操作每個(gè)單元都需要.間接計(jì)算方法主要是通過建立各種規(guī)則化邊界積分方程,來間接計(jì)算CPV和HFP積分.包括直接變量規(guī)則化邊界積分方程和間接變量規(guī)則化邊界積分方程[7-9,15-20].
不同于已有的工作,本文提出了一個(gè)新的計(jì)算強(qiáng)奇異積分的方法,稱為簡單解法.它利用簡單解過程,獲得強(qiáng)奇異積分值,無需直接計(jì)算積分,也無需建立規(guī)則化邊界積分方程,因此方法具有理論簡單、計(jì)算效率高、結(jié)果準(zhǔn)確等優(yōu)點(diǎn).此外,方法能夠計(jì)算任何邊界通量?u/?xi(i=1,2).該方法容易跟隨和拓展到其他如彈性力學(xué)問題、Stokes問題、Helmholtz問題等.
本文假定Ω是R2中的一個(gè)有界區(qū)域,Ωc是它的補(bǔ)域,Γ=?Ω是它們的邊界.n(x)是區(qū)域Ω的邊界Γ在x點(diǎn)處的單位外法向量.位勢問題的基本解是
(1)


(2)
Ω內(nèi)的位勢可表示為[17]

(3)
(4)
對方程(4),令y→Γ,可獲得如下邊界積分方程
(5)

(6)
為了能夠更清楚地表述問題的本質(zhì),采用方程(6)的離散化形式.將邊界Γ分割成N段(Segment),并假定y在第l段上,即y∈Γl,并記為yl,于是有
(7)
式(7)右端第一項(xiàng)中的積分均為正常積分,因此問題的關(guān)鍵在于如何估計(jì)式(7)右端第二項(xiàng)積分.若記
方程(7)可寫為
φ(yl)A,yl∈Γl
(8)
當(dāng)常元插值(常單元)被采用時(shí),方程(8)可寫為
yl∈Γl
(9)
這里φm表示密度函數(shù)φ在第m個(gè)節(jié)點(diǎn)上的值.為了確定(7)式中的A,提出一種簡單、高效、準(zhǔn)確的簡單解法:
(1)對于有限域Ω,引入如下簡單解
u(x)=x1+x2+1
(10)
(2)對于無限域Ωc,引入如下簡單解
(11)
這里(a,b)是Ωc外的任一點(diǎn).
注:原則上,簡單解可以任意選取.對于內(nèi)邊值問題,簡單解只需滿足Laplace方程;對于外邊值問題,簡單解除滿足Laplace方程外,還需滿足無窮遠(yuǎn)邊界條件.
方程(4)的離散化形式
l=1,…,N
(12)
這里φm是邊界密度函數(shù)φ(x) 在第m個(gè)單元Γm(m=1,…,N)上的節(jié)點(diǎn)處的值.
對于有限(無限)域問題,將簡單解式(10)(或式(11))代入方程(12)中,可求得密度函數(shù)在N個(gè)節(jié)點(diǎn)處的值φm,m=1,…,N.然后,將所求得的密度函數(shù)解φm及簡單解式(10)(或式(11))代入方程(9)中,即可求得A.
以上方程經(jīng)過適當(dāng)?shù)慕M合變形后,可適用于任何邊界條件的邊值問題.
利用三個(gè)數(shù)值算例,以驗(yàn)證本文算法的精確性、有效性及收斂性.所有算例中,邊界離散單元的邊界量采用常元插值,單元幾何邊界采用精確單元L2范數(shù)下的相對誤差定義為
Relative error =
(13)

算例1首先考察一個(gè)具有Dirichlet不連續(xù)邊界條件的的穩(wěn)定溫度場問題.計(jì)算域是[0,1]×[0,1],邊界條件是
(14)
問題的解析解是
u(x1,x2)=
(15)
其中
(16)
計(jì)算時(shí),邊界單元數(shù)是80個(gè).圖1(a)和圖1(b)分別描述了解析法和本文數(shù)值法獲得的場位勢的等值圖,可觀察到,本文數(shù)值解和解析解比較吻合.為了進(jìn)一步考察本文方法的準(zhǔn)確性和有效性,選擇1 600個(gè)計(jì)算點(diǎn)均布在區(qū)域0.15≤x,y≤0.85內(nèi),圖2給出了這些計(jì)算點(diǎn)上的場溫度數(shù)值解的絕對誤差曲面.從中可以看出,盡管具有很少的邊界單元數(shù),本文獲得的數(shù)值結(jié)果仍然相當(dāng)準(zhǔn)確,表明本文給出的方法能夠有效地求解具有不連續(xù)Dirichlet邊界條件的邊值問題.意味著,本文計(jì)算強(qiáng)奇異積分對角元的方法是非常有效的.

(a) 解析解 (b) 數(shù)值解圖1 場位勢解Fig.1 The field potential solution

圖2 邊界通量解的收斂曲線Fig.2 Convergence curves of boundary flux solutions


(17)

圖3 計(jì)算域Fig.3 Problem sketch
計(jì)算時(shí),外邊界Γ1離散為36單元,內(nèi)邊界離散為24個(gè)單元.表1給出了沿x軸一些域內(nèi)點(diǎn)處溫度的數(shù)值結(jié)果;圖4(a)和(b)分部描述了邊界Γ1上的法向通量解和邊界Γ2的溫度解的相對誤差.從中可看出,即使使用很少的單元數(shù),數(shù)值解已相當(dāng)?shù)販?zhǔn)確.
表1 內(nèi)點(diǎn)溫度解
Tab.1 Temperature solutions on internal points

內(nèi)點(diǎn)精確解數(shù)值解相對誤差(6.0, 0.0)0.1400000×1030.1399699×1032.149503×10-4(7.0, 0.0)0.1840000×1030.1839418×1033.163785×10-4(8.0, 0.0)0.2340000×1030.2339130×1033.719932×10-4(9.0, 0.0)0.2900000×1030.2898743×1034.333464×10-4

(a) 通量q

(b) 溫度u圖4 邊界量數(shù)值解的相對誤差Fig.4 Ralative error of boundary temperature uand flux q
算例3考慮不規(guī)則區(qū)域上,具有混合邊界條件的溫度場問題,如圖5所示.邊界Γ由兩部分構(gòu)成,即Γ=Γ1∪Γ2,其中

圖5 計(jì)算域Fig.5 Problem sketch



問題的解析解是如下振蕩調(diào)和函數(shù)
u(x1,x2)=ex1+3cosx2
(18)
圖6 給出了邊界Γ上的位勢函數(shù)u的圖像,從圖中可觀察到,u的最大值大于300,最小值小于-20,表明位勢邊界函數(shù)u具有很強(qiáng)的振蕩性.

圖6 邊界Γ上的位勢函數(shù)Fig.6 Potential boundary function on Γ
20個(gè)域內(nèi)計(jì)算點(diǎn)均布在中心在(0.5,0.5),半徑是0.4的圓周上,如圖5所示.圖7描述了這些計(jì)算點(diǎn)上的溫度數(shù)值解u的收斂曲線,圖8給出了邊界Γ1上的通量和Γ2上的溫度數(shù)值解的收斂曲線.從中可看出,不論是內(nèi)點(diǎn)還是邊界點(diǎn),隨著邊界單元數(shù)的增加,相對誤差迅速減小,可見本文方法在求解具有振蕩解析的問題時(shí)仍具有很好的收斂特征.表明本文計(jì)算強(qiáng)奇異對角元的方法是非常有效的.

圖7 內(nèi)點(diǎn)溫度數(shù)值解的收斂曲線Fig.7 Convergence curves of temperature solutions at internal points

圖8 邊界溫度u和通量q的收斂曲線Fig.8 Convergence curves of boundary temperature uand flux q
本文提出了計(jì)算強(qiáng)奇異積分的簡單解法,適用于任何邊界通量的計(jì)算,具有數(shù)學(xué)理論簡單,程序設(shè)計(jì)容易,計(jì)算準(zhǔn)確性高,容易跟隨和推廣等優(yōu)點(diǎn).數(shù)值算例驗(yàn)證了方法的有效性,表明本文方法能夠有效準(zhǔn)確地計(jì)算強(qiáng)奇異積分的對角元.