邊界元法的應用

廉蒙蒙

（安陽工學院機械工程學院，河南安陽455000）

0 引言

邊界元法是用于求解工程與科學問題的常用數值分析方法之一，它的基本思想是將經典問題的基本解融入到邊界積分方程中，把微分方程的邊值問題歸化成邊界上的積分方程，并以其為控制方程，在借鑒有限元離散的基礎上，利用各種離散化技術，通過將邊界進行插值離散后，邊界積分方程離散為代數方程組,再用數值方法求解代數方程組,從而得到原問題邊界積分方程的解[1-2]。邊界元法具有邊界積分方程非常深厚的數學基礎,是繼有限元法之后興起的一種有效的數值方法，這種方法不僅可以用來分析彈性力學等固體力學的問題,同時也可應用于流體力學、聲場,以及電磁場等其它的研究領域[3]。

早在1963年M.A.Jaswon[4]首次通過有勢的理論對Laplace方程變換，建立了邊界積分方程的數值方法，即誕生了間接邊界元法，并求解了位勢問題。同年G.T.Symm也做了相關的工作[5]。Rizzo在1967年首次提出了直接邊界元法，并求解了二維彈性靜力學問題，這也是邊界元法文獻中第一篇關于直接邊界元法的文獻[6]。關于邊界積分方程法的專著是T.A.Cruse和F.J.Rizzo[7]于1975年出版的，并正式的命名此法為邊界積分方程方法，這也是邊界元中的第一本專著。C.A.Brebbia[8]于1986年編著了《工程師用的邊界單元法》，全面詳細地介紹了邊界元法的原理和應用，并討論了位勢問題和二維彈性靜力學力學問題在邊界元方法中的應用，并用FORTRAN語言編輯的程序來數值模擬相關問題，這標志著邊界元方法已進入了系統性的研究。

1 邊界積分方程

邊界元法中運用Betti功能互等定理，得到著名的Somigliana等式，即區域內點Y位移的積分表達式[9]

當Y→Γ時，得到彈性靜力學的位移邊界積分方程

式中 ∫Γ(? ? ?)dΓ為沿邊界積分的Cauchy主值積分，x和y分別表示邊界上的場點和源點，如圖1所示。Cij(y)稱為位移奇異系數，它依賴于邊界Γ上y點附近的幾何形狀與彈性常數，對于光滑的邊界點處，Cij=δij2。但是在實際的計算中，并不是直接計算這一系數，而是通過剛體位移求得。uj(x)和tj(x)分別表示邊界Γ上的位移分量和面力分量；Ui?j(x,y)，Ti?j(x,y)為邊界積分方程的核函數，亦稱為基本解，它們分別表示在無限彈性平面的一點y(源點)的i方向作用單位集中力時，在邊界上ui任意一點x(場點)的j方向上產生的位移和面力。對于平面應變問題它們的具體形式為[10]

式中，令xi和yi分別為場點和源點的坐標,則源點y到場點x的矢徑r→定義為

采用線性單元離散邊界積分方程后，可以寫成

把各節點的位移與面力分別記作單元位移矢量{u(e)}和單元面力矢量{t(e)}

2數值計算

例：受均勻壓力P=1的狹長薄板，如圖2所示，狹長板沿x軸的長度為L=2m，板的厚度h在1～10-10m之間變化，板的彈性剪切模量μ=8.0×1010Pa，泊松比υ=0.2，當y=0時，邊界條件ux=uy=0。

圖2 受均勻壓力的狹長板

圖3 A點的應力

圖4 A點的應力σy

圖3和4的數據表明，本文所計算的結果即使是在h=10-10m時，依然與精確解非常的吻合。表1、表2的數據表明，本文即使沒有進行近奇異積分的處理，但是所計算的結果依然和精確解十分吻合。

3 結論

本文應用邊界元法所計算的結果即使沒有處理近奇異積分，依然比較精確，究其原因是Mathe?matica軟件數值精度比較高，其次是智能積分，通常會采用自適應算法計算積分的近似值，對積分區間進行遞歸分割，直到達到指定的準確度為止。以往近奇異積分的處理即坐標變換是為了將變化比較劇烈的曲線變得比較平緩，這樣用Gauss積分法計算近奇異積分式子才能達到很高的精度。由此可見，計算近奇異積分，對于高級的程序語言來說，其智能的積分計算功能已經完全能夠勝任近奇異積分計算這項工作，不需要再做任何的處理。

表1 直線x=1.8上內點應力σx的計算結果(h =10-9m)