自然邊界元與有限元求解平面彈性問題的耦合法

1、自然邊界元與有限元求解平面彈性問題的耦合法臧彤,趙慧明,楊敏(中國(guó)礦業(yè)大學(xué)力學(xué)與建筑工程學(xué)院,江蘇 徐州 221008摘要:為了更充分地利用有限元與自然邊界元各自的優(yōu)點(diǎn),并盡量減少由于方法的局限性造成的在計(jì)算量及計(jì)算精度上的不足。本文引入了有限元與自然邊界元耦合的方法,在文中簡(jiǎn)單介紹了自然邊界元法及其與有限元的耦合的原理,通過設(shè)置含有重疊區(qū)域的圓形人工邊界,實(shí)現(xiàn)自然邊界元法與有限元法的耦合,并把該方法應(yīng)用到無界區(qū)域上的實(shí)際算例中,從計(jì)算結(jié)果的比較中可以看出自然邊界元與有限元耦合算法的收斂速度更快,充分體現(xiàn)了耦合法在解決無界區(qū)域問題上的優(yōu)越性。關(guān)鍵詞:自然邊界元法;有限元法;無界區(qū)域;耦合法;D

2、-N迭代中圖分類號(hào):TB112A coupling method of natural boundary element and finiteelement for the planar elastic problemZang Tong, Zhao Huiming, Yang Min(School of Mechanices and Civil Engineering, China University of Mining and Technology,JiangSu XuZhou 221008Abstract: In order to make better use of finite el

3、ement and boundary elements the merits, and to minimize the deficiencies in calculating the amount and accuracy on the deficiencies caused by method limitations, this article introduced the natural boundary element method and its coupling with the finite element process, and described simply the nat

4、ural boundary element method and its coupling with the finite element process in the text. This paper gave a briefing on the natural boundary element iterative realization of D-N by setting a circular artificial boundary containing the overlap region, and applied the method to the real unbounded exa

5、mple. By the comparing of results from the different calculation methods, the precision and algorithm converges of results by using coupling method were more good. The analysis fully reflected the superiority of coupling method to solve unbounded regional problems. Keywords:NBEM; FEM; unbounded regi

6、onal; direct coupling method; D-N iteration0引言自然邊界元法與有限元法自創(chuàng)立以來,在眾多領(lǐng)域都取得了矚目的成就。但它們也都存在著一定的不足,有限元方法在處理無界問題時(shí)有著本質(zhì)的困難,而自然邊界元法則由于一般區(qū)域的格林函數(shù)難以求得,使得它的應(yīng)用范圍較窄,但二者基于相同的變分原理,這為兩種方法的耦合奠定了堅(jiān)實(shí)的理論基礎(chǔ)。自然邊界元方法,是由Green函數(shù)和Green公式出發(fā),將微分方程邊值問題歸化為邊界上強(qiáng)奇異積分方程(或稱為超奇異積分方程的一種數(shù)值計(jì)算方法,在一定程度上它保持了原邊值問題的一些有用的性質(zhì),保證了解的唯一性和穩(wěn)定性【1】。另外從數(shù)值計(jì)算的

7、角度,與一般邊界元方法相比,自然邊界元方法使剛度矩陣的計(jì)算量也大為減少。耦合法既吸取了有限元能適應(yīng)較任意的有限區(qū)域的優(yōu)點(diǎn),又克服了自然邊界歸化對(duì)區(qū)域一些要求的限制,其總體剛度矩陣正是自然邊界元法及有限元法得到的剛度矩陣之和【2】。這與其他類型的邊界元與有限元的耦合法相比要簡(jiǎn)單的多。本文借助于區(qū)域分解算法的思想,并基于圓形區(qū)域人工邊界的自然邊界歸化理論和有限元基本理論,運(yùn)用自然邊界元與有限元的耦合法實(shí)現(xiàn)對(duì)無界區(qū)域上含有重疊區(qū)域的平面彈性問題的分析及求解。本文將通過有重疊區(qū)域迭代算法實(shí)現(xiàn)自然邊界元與有限元的耦合,并將該方法應(yīng)用于求解帶方孔的無界彈性區(qū)域問題。1 耦合法在無界區(qū)域問題中的重疊型迭代實(shí)

8、現(xiàn)原理自然邊界元與有限元耦合法在求解無界區(qū)域問題時(shí)尤其顯示其優(yōu)越性。由于有限元法往往難以處理無界區(qū)域問題。故在工程計(jì)算中常引入人工邊界,以割去區(qū)域的無界部分而僅在剩下的有界區(qū)域中求解【3】。這樣為滿足足夠的精度,必須取相當(dāng)大的求解區(qū)域,而計(jì)算一個(gè)相當(dāng)大的區(qū)域上的邊值問題的數(shù)值解依然是很困難的。自然邊界元與有限元的耦合法可克服這一困難。通常取圓周為人工邊界,即將區(qū)域分為圓周外部區(qū)域及剩下的有界區(qū)域,對(duì)圓外區(qū)域應(yīng)用自然邊界歸化。于是原邊值問題化為有界區(qū)域上的等價(jià)變分問題,然后用自然邊界元與有限元的耦合數(shù)值求解之【4】。 圖1 區(qū)域分解圖如圖,設(shè)為包圍原點(diǎn)的充分光滑的有界閉曲線,作其外區(qū)域,為包圍的

9、,以原點(diǎn)為圓心,以R 為半徑的邊界元區(qū)域2的邊界。為包圍的,以原點(diǎn)為圓心,以r 為半徑的有限元區(qū)域1的邊界。本文在具體實(shí)現(xiàn)自然邊界元與有限元的耦合法時(shí)采用了迭代法來實(shí)現(xiàn),下面介紹一下無界平面彈性問題的自然邊界元與有限元的耦合法的有重疊區(qū)域迭代原理?？紤]有光滑邊界的平面區(qū)域上的平面彈性方程組的第一邊值問題=+上內(nèi),0(0GG GG graddiv (1 如上圖所示, 1為有界區(qū)域,2為無界的圓外區(qū)域【5】。1區(qū)域用有限元方法來進(jìn)行求解,2用自然邊界元方法來求解,含重疊型區(qū)域問題的有重疊區(qū)域迭代法的實(shí)現(xiàn)過程如下:(1 建立模型并根據(jù)實(shí)際情況或有限元計(jì)算結(jié)果將其劃分為1有限區(qū)域、2無限區(qū)域,區(qū)域1與

10、2有重疊。(2 給出2區(qū)域邊界上的初始位移0u G(3 利用Poisson 積分公式(=,*,R u R u P PP P r u r u r rr rrr (R r > (3 計(jì)算出有限元邊界上的位移(,r u r和(,r u 【6】。(4 將計(jì)算得到的(,r u r 和(,r u 作為1區(qū)域的邊界條件,利用有限元公式,解方程組P u K GG =其中K 為有限元?jiǎng)偠染仃囉?jì)算得到位移u G ,P G為1有限區(qū)域初始應(yīng)力條件。(5 將u G 與0u G比較,如果其差的絕對(duì)值滿足精度要求,則結(jié)束計(jì)算,否則進(jìn)行第六步(6 再將上步u K作為2區(qū)域的位移邊界條件0u G ,代入第二步計(jì)算,從而

11、得到第二次循環(huán)下的u K ,再進(jìn)行二者的比較,如此反復(fù),直到誤差達(dá)到要求。2 數(shù)值算例2.1 算例例:考慮研究帶方孔的無界彈性區(qū)域,方孔是正方形,邊長(zhǎng)為4m ,方孔周邊施加1000m N /的均布載荷,且正方形的形心正好在原點(diǎn)處。設(shè)材料的彈性模量E=40GPa ,泊松比為0.3。模型如下圖所示:我們運(yùn)用上述的有重疊區(qū)域迭代算法來實(shí)現(xiàn)該此問題的耦合。 圖(2 截面劃分圖選取r 為定值240m ,則選取不同的R 值時(shí),自然邊界元與有限元耦合法與有限元法在點(diǎn)A (如上圖即(2,2點(diǎn)的X 方向位移和Y 方向位移值如下表所示。并對(duì)兩種方法得到的位移值進(jìn)行比對(duì)分析,變化曲線如下圖: 圖(3A(2,2點(diǎn)的x

12、方向位移圖 圖(4A(2,2點(diǎn)的y方向位移根據(jù)計(jì)算結(jié)果畫出R=20m時(shí)的主應(yīng)力等值線,如下圖所示: 應(yīng)力等值線 圖(6 ANSYS 下應(yīng)力圖2.2 松弛因子的影響為了加快迭代的速度,引入松弛因子。松弛因子即:在每一步迭代時(shí),計(jì)算出kU ,對(duì)其進(jìn)行一次松弛,k k k U U U1('1('+=【7】;的取值范圍為01,下面我們就此問題通過實(shí)驗(yàn)給出表3表3 不同松弛因子的收斂速率比較松弛因子 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 迭代次數(shù) 12 9 7 6 8 11 15 23 35注:松弛因子取0.4時(shí),收斂速度最快。2.3 結(jié)果分析在應(yīng)用有限

13、元法解決上述算例時(shí),在半徑比較小的時(shí)候,收斂結(jié)果不是很好。 有限元的收斂速度與劃分網(wǎng)格的疏密有關(guān),為了控制計(jì)算量,計(jì)算時(shí)網(wǎng)格劃分的較稀疏,在一定程度上影響了計(jì)算結(jié)果的精度。3 結(jié)論通過前面的研究得到如下結(jié)論(1用自然邊界元法解決區(qū)域邊界問題,劃分的網(wǎng)格較少,減小了計(jì)算量,并且保證了精度。在半徑比較小的時(shí)候,耦合法就得到了比較精確的結(jié)果,并且隨著重合區(qū)域?qū)挾鹊脑黾?位移收斂情況更好。(2松弛因子的選取對(duì)迭代收斂速度的影響很大。當(dāng)松弛因子取0.4時(shí),迭代收斂速度最快。(3另外,本文只是針對(duì)有重疊區(qū)域進(jìn)行的耦合,還可對(duì)耦合區(qū)域?yàn)闊o重疊區(qū)域的情況進(jìn)行進(jìn)一步研究。 參考文獻(xiàn) (References1 余德浩.自然邊界元方法的數(shù)學(xué)理論M.北京.科學(xué)出版社,19932 余德浩.