李威，李諍，柴應(yīng)彬，韓旭

（1.華中科技大學船舶與海洋工程學院，武漢430074；2.湖南大學機械與運載工程學院，長沙410082）

基于邊光滑有限元法的結(jié)構(gòu)振動與穩(wěn)定性研究

李威1，李諍1，柴應(yīng)彬1，韓旭2

（1.華中科技大學船舶與海洋工程學院，武漢430074；2.湖南大學機械與運載工程學院，長沙410082）

用一種新型的數(shù)值方法—邊光滑有限元方法（edge-based smoothed finite element method,ES-FEM）對梁、板等簡單結(jié)構(gòu)的振動與穩(wěn)定性問題進行分析與研究。光滑有限元方法在計算系統(tǒng)總體剛度矩陣時在構(gòu)造的光滑域內(nèi)對應(yīng)變進行光滑操作，這樣得到的剛度矩陣改善了有限元方法過剛的數(shù)值缺陷。對梁和板結(jié)構(gòu)在一定邊界條件下的振動與穩(wěn)定性問題的分析結(jié)果表明：與其它有限元方法相比，光滑有限元方法提高了數(shù)值計算結(jié)果的精確度和收斂性。

振動與波；有限元法；邊光滑有限元法；剛度矩陣；穩(wěn)定性

Kay words:vibration and wave;FEM;ES-FEM;stiffness matrix;stability

梁和板結(jié)構(gòu)是機械、車輛、土木、船舶等工業(yè)領(lǐng)域最基本的構(gòu)件，它們的振動與穩(wěn)定性關(guān)系到設(shè)備及整個系統(tǒng)的工作性能和壽命，嚴重時會使系統(tǒng)功能失效或者造成破壞事故，因此，對結(jié)構(gòu)的振動和穩(wěn)定性研究有重要的現(xiàn)實意義，而有限元方法作為重要的數(shù)值方法被廣泛應(yīng)用，其中三角形單元以其方便的前處理和對復雜幾何形狀的適應(yīng)性，深受研究者的喜愛。然而，有限元三角形單元的發(fā)展深受精度和穩(wěn)定性低的制約。應(yīng)變光滑技術(shù)是CHEN[1]在節(jié)點積分無網(wǎng)格方法中為得到穩(wěn)定解提出來的，很快被應(yīng)用到自然單元方法。LIU[2―4]針對有限元方法(FEM)“過剛”的問題，將標準有限元方法和應(yīng)力光滑技術(shù)結(jié)合起來提出了基于單元光滑有限元方法(SFEM)和基于節(jié)點光滑有限元方法(NS-FEM).與FEM“過剛”的特性正相反，NS-FEM表現(xiàn)出“過柔”的問題，而過柔的模型將會導致結(jié)果的不穩(wěn)定性。針對此問題，LIU[5]提出一種基于邊光滑有限元方法用于計算平面問題，展示了較好的計算性能。

本文基于邊光滑有限元方法(ES-FEM)對梁、板等簡單結(jié)構(gòu)進行振動和穩(wěn)定性分析,先將問題域離散為三角形單元，基于三角形的邊進一步構(gòu)造光滑域，在光滑域內(nèi)運用應(yīng)變光滑技術(shù)調(diào)整系統(tǒng)剛度，改善有限元方法過剛的模型缺陷，提高計算精度。

1 邊光滑有限元法基本公式

1.1 光滑域的構(gòu)造

在邊光滑有限元法中，問題域Ω首先被離散成N個三角形單元，將每條邊的端點和這條邊相鄰的兩個單元的中心相連接，這樣就在三角形單元的基礎(chǔ)上形成了Nk個基于邊的光滑域，如圖1所示。

圖1 三角形單元和基于邊的光滑域

1.2 靜力分析

對于平面應(yīng)力問題，ES-FEM的靜力平衡方程[5]為

1.3 動力分析

平面問題的動力平衡方程可表示為

把d=d0exp(iwt)代入上式，不計阻尼和外力，可得到特征方程

1.4 Reissner-mindlin板問題

設(shè)板離散成若干個三角形單元，每個單元有三個節(jié)點，每個節(jié)點有三個位移：每個節(jié)點的位移可表示為

自由振動分析的標準Galerkin弱形式[6]可表示為

對于板的小撓度屈曲，標準Galerkin弱形式[7]可表示為

由標準Galerkin弱形式可得板靜力問題和屈曲問題的ES-FEM公式，導出特征值方程分別為

式中NS為光滑域的數(shù)目，為光滑域的面積，Ne為初始單元的數(shù)目，τ為板內(nèi)初始應(yīng)力矩陣，并且有

式中[8]k=5/6為剪切修正系數(shù)，α=0.05，h為三角形單元最長邊的邊長。對和有

式中Nj為第j個節(jié)點的形狀函數(shù)。

2 數(shù)值算例

2.1 懸臂梁自由振動分析

考慮一懸臂梁模型，長L=120 mm,高H=12 mm,厚度t=1.0 mm，楊氏模量E=2.1×1011N/m2，泊松比υ=0.3，密度ρ=8×10-6kg/mm3。用單元尺寸為0.001m個節(jié)點的標準有限元方法結(jié)果作為參考，用Matlab軟件編程實現(xiàn)FEM和ES-FEM，在單元尺寸為0.006 m和0.004 m的情況下分別對算例進行求解，前5階固有頻率如表1所示。

表1 懸臂梁自由振動前5階固有頻率

從表1可以很明顯的觀察到，當用同樣的節(jié)點數(shù)目時，ES-FEM可以得到比FEM更精確的結(jié)果，收斂的更快。固有頻率可以用來衡量結(jié)構(gòu)的剛度，因此與FEM相比，用ES-FEM建立的模型更接近結(jié)構(gòu)的精確剛度。用ES-FEM方法求懸臂梁自由振動前八階模態(tài)振型如圖2所示

2.2 Reissner-mindlin板自由振動分析

有一兩端固支，兩端自由的方形板，邊長L=0.5 m，板厚t=0.002 5 m，材料楊氏模量E=2.0×1011N/m2，泊松比υ=0.3，質(zhì)量密度ρ=8 000 kg/m3，采用2×20×20個三角形單元對板進行自由振動分析，用無量綱量表示結(jié)構(gòu)的固有頻率，這里D=Et3/(12(1-υ2))為板的彎曲剛度。求解結(jié)果和精確解[9]如表2所示。

圖2 懸臂梁前8階模態(tài)振型

表2 方板自由振動前4階無量綱頻率

同精確解相比，在同樣的節(jié)點數(shù)和單元劃分條件下，ES-FEM方法得到的解比FEM更精確。用ESFEM方法求方板自由振動的前4階模態(tài)振型如圖3所示。

2.3 Reissner-Mindlin板穩(wěn)定性分析

考慮一方形板，四邊固支，單軸向受力，邊長L=1 m，板厚t=0.01 m，材料楊氏模量E=2.0×1011N/m2，泊松比υ=0.3，定義屈曲載荷系數(shù)K=λcrL2/(π2D)，D為前例定義的彎曲剛度。分別用尺寸為0.1 m和0.05 m的三角形單元對板進行穩(wěn)定性分析，解析解[10]為K=10.07，數(shù)值方法計算結(jié)果如表3所示。

表3 方形板屈曲失穩(wěn)載荷系數(shù)

用ES-FEM求得四邊固支方形板軸向失穩(wěn)模態(tài)如圖4所示。

由FEM和ES-FEM的計算結(jié)果和解析解對照，可以明顯的看出，在同樣的單元劃分條件下，與FEM相比，ES-FEM更接近精確解，在增加單元時收斂更快。

3 結(jié)語

本文采用一種新的邊光滑有限元方法對梁、板等結(jié)構(gòu)進行穩(wěn)定性分析，在解析解的參照下將結(jié)果與標準有限元方法進行對比，結(jié)果表明：邊光滑有限元方法有效地改善了標準有限元方法過剛的模型缺陷，在同樣的單元條件下結(jié)果更為精確，收斂性更快更好，是一種比較理想的數(shù)值計算方法。

[1]CHEN J S,WU C T,YOON S,et al.A stabilized conformingnodalintegrationforGalerkinmeshfree method[J].International Journal for Numerical Methods in Engineering,2000,50:435-466.

[2]LIU G R,Nguyen T T,Dai K Y,et al.Theoretical aspects of the smoothed finite element method(SFEM)[J].InternationalJournalforNumericalMethodsin Engineering.2007,71(8):902-930.

[3]LIU G R,DAI K Y,Nguyen T T,et al.A smoothed finite elementmethodformechanicsproblems[J].Computational Mechanics,2007,39(6):859-877.

[4]LIU G R,A generalized gradient smoothing technique and smoothed bilinear form for Galerkin formulation of a wide class of computational methods[J].International Journal of Computational Methods.2008,5(2):199-236.

[5]LIU G R,DAI K Y,Nguyen T T.An edge-based smoothed finite element method(ES-FEM)for static,free and forced vibration analysis[J].Journal of Sound and Vibration, 2009,320:1100-1130.

[6]LIU G R,Quek S S.The finite element method[M]:A practical course.Butterworth Heinemann,Oxford,2002.

[7]Liew K M,CHEN X L.Buckling of rectangular Mindlin plates subjected to partial in-plane edge loads using the radial point interpolation method[J].International Journal of Solids and Structures.2004,41:1677-1695.

[8]Lyly M,Stenberg R,Vihinen T.A stablebilinear element for the Reissner-Mindlin plate model[J].ComputerMethodsinAppliedMechanicsand Engineering,1993,110:343-357.

[9]Leissa A W.Vibration of plates[M].Washington,DC, NASA,SP-160,1969.

[10]Timoshenko S P,Gere J M.Theory of elastic stability [M].3rd edition.New York,McGraw-Hill,1970.

[11]曾攀.有限元分析及應(yīng)用[M].北京：清華大學出版社，2004.

Analysis of Vibration and Stability of Structures Using Edge-based Smoothed Finite Element Method

LIWei1,LIZheng1,CHAI Ying-bin1,HANXu2

(1.School of NavalArchitecture and Ocean Engineering,
Huazhong University of Science and Technology,Wuhan 430074,China; 2.College of Mechanical and Vehicle Engineering,Hunan University,Changsha 410082,China)

A novel numerical method,edge-based smoothed finite element method(ES-FEM),is used to analyze vibration and stability problems of beam and plate structures.In the ES-FEM,the system global stiffness matrix is computed using strain-smoothed method over the smoothing domains connected to the edges of the elements.As a result,the smoothing operation can reduce the over-stiffening behavior of the conventional finite element method.Numerical examples of beam and plate structures under special boundary conditions are presented to demonstrate that the ES-FEM can achieve higher accuracy and convergence rate compared with the conventional FEM for structures.

O241.82

ADOI編碼：10.3969/j.issn.1006-1335.2014.02.001

1006-1355(2014)02-0001-04

2013-01-09

國家杰出青年科學基金項目（10725208）

李威（1975-），男，湖北武漢人，副教授，目前從事結(jié)構(gòu)振動噪聲分析與試驗測試研究。

E-mail:hustliw@mail.hust.edu.cn