羅 瑪，張大義，朱海雄，馬艷紅,3，洪 杰,3

（1.中國航發(fā)湖南動力機械研究所，湖南株洲412002；2.北京航空航天大學能源與動力工程學院，北京100191；3.先進航空發(fā)動機協(xié)同創(chuàng)新中心，北京100191）

基于試驗數(shù)據(jù)識別的金屬橡膠有限元仿真方法

羅 瑪1，張大義2，朱海雄2，馬艷紅2,3，洪 杰2,3

（1.中國航發(fā)湖南動力機械研究所，湖南株洲412002；2.北京航空航天大學能源與動力工程學院，北京100191；3.先進航空發(fā)動機協(xié)同創(chuàng)新中心，北京100191）

金屬橡膠工程應用分析常需要1種可以與商用有限元軟件相結合，適用于復雜結構設計的快捷數(shù)值仿真方法。針對金屬橡膠的有限元仿真問題，基于A N SY S軟件，發(fā)展了基于試驗數(shù)據(jù)識別的橫向各向同性材料參數(shù)確定方法，并采用循環(huán)加載、模型更新的方法，實現(xiàn)了金屬橡膠材料受壓過程非線性力學特性的仿真計算。結果表明：該方法對立方體件在無約束狀態(tài)下受壓過程進行模擬時，相對誤差在10%之內(nèi)；安裝約束使其剛度增大，仿真與試驗結果相比規(guī)律吻合，剛度的誤差在30%之內(nèi)；環(huán)形金屬橡膠件在徑向受壓時，仿真與試驗測得剛度的相對誤差在25%之內(nèi)。

金屬橡膠；力學性能；有限元；仿真計算；復雜受力狀態(tài)

0 引言

金屬橡膠是1種性能優(yōu)良的阻尼材料，在航空、航天、船舶等領域的振動抑制方面有著廣泛應用[1-2]。其力學性能的試驗標定與理論研究一直是國內(nèi)外學術界和工程界備受關注的研究熱點。

在試驗研究方面，Gadot[3]、Li[4]、馬艷紅[5]、洪杰[6]、Liu[7]、Zarzour[8]、Zhang[9-10]等分別研究得到了金屬橡膠的結構參數(shù)、材料參數(shù)和工藝參數(shù)對其剛度特性和阻尼性能的影響規(guī)律，為金屬橡膠的工程應用提供了充分的數(shù)據(jù)支撐。

在力學模型構建方面，李宇燕等[11]提出了金屬橡膠材料的非線性干摩擦結構單元力學模型；于鋒禮等[12]基于黏彈性阻尼機理提出了黏彈性本構關系模型；張玲凌等[13]提出了將金屬橡膠的遲滯回線分解為彈性力和阻尼力的等效法，從而建立了恢復力分解模型；朱彬等[14]提出了金屬橡膠的螺旋卷模型；Ma等[15]構建了可以同時對金屬橡膠彈性性能和耗能特性進行定量描述的“傾斜螺旋微元體+摩擦角錐”的細觀結構單元體；王亮[16]通過不可壓縮、各向同性和力學性能參數(shù)假設，建立了金屬橡膠密封系統(tǒng)的3維實體接觸模型。

上述研究工作雖為金屬橡膠的力學性能預測奠定了良好的理論基礎，但在工程實際中，常常需要1種適用于復雜結構設計的快捷數(shù)值仿真方法，并且可以與商用有限元軟件相結合，從而滿足工程人員的分析需求，本文針對此需求開展研究工作。

1 金屬橡膠材料參數(shù)及確定方法

1.1 關鍵參數(shù)確定

金屬橡膠的特殊成型工藝，導致其材料力學性能參數(shù)具有橫向各向同性的特征[17-18]，即成型方向的材料參數(shù)特性與非成型方向不同，各非成型方向上則相同。此外，其成型方向和非成型方向上的力學性能參數(shù)在大應變范圍內(nèi)體現(xiàn)出明顯的非線性，彈性模量、泊松比和剪切模量均表現(xiàn)出明顯的應變相關性[9，19]。

為明確定義各參數(shù)，下標1表示成型方向，2表示非成型方向，12表示由成型方向和非成型方向組成的材料參數(shù)面，22表示由2個非成型方向組成的材料參數(shù)面。對于金屬橡膠，按各向異性材料模型給定各參數(shù)項時，應包括 E1、E2、v12、v22、G12、G22。

由于2個非成型方向組成的材料參數(shù)面為各向同性面，因此面內(nèi)的彈性模量、泊松比和剪切模量符合廣義胡克定律，有如下關系式

因此，金屬橡膠的獨立材料參數(shù)僅有E1、E2、v12、v22、G12，確定這5個參數(shù)，就可以建立仿真計算所需的模型。

1.2 參數(shù)確定方法

采用準靜態(tài)壓縮試驗測試金屬橡膠標準試件在各不同材料方向上的力學性能，試驗采用斜坡位移控制方法，每秒應變加載0.001。以相對密度為0.21，幾何尺寸為20 mm×20 mm×20 mm的金屬橡膠構件為例，給出試驗設備及試件，如圖1所示。根據(jù)試驗測試結果，不同方向上彈性模量隨應變的變化曲線如圖2所示。采用最小二乘法擬合曲線，得到彈性模量隨應變的函數(shù)關系式

圖1 準靜態(tài)試驗設備及正方體試件

圖2 各彈性模量隨應變變化曲線

在不同方向上受載時，測得各泊松比隨應變的變化曲線，如圖3所示。經(jīng)最小二乘法擬合得到其分別隨應變變化的關系式

在成型方向上受載時，剪切模量隨應變的變化關系曲線如圖4所示。曲線擬合后得到其隨應變變化的關系式

圖3 各泊松比隨應變變化曲線

圖4 剪切模量G12隨應變變化曲線

在確定金屬橡膠材料力學性能參數(shù)時，按照橫向各向同性材料模型確定，即給定以上各參數(shù)表達式，并將表達式中的應變與材料模型在不同方向上的應變一一對應。

2 有限元仿真模型與計算方法

2.1 仿真計算模型

基于ANSYS有限元軟件，采用solid185單元，分別針對立方體件(20 mm×20 mm×20 mm)和環(huán)形件(內(nèi)/外徑 100 mm/116 mm，厚度 10 mm)，建立金屬橡膠的仿真計算模型，如圖5所示。材料參數(shù)按各向異性材料模型確定，各方向力學性能參數(shù)按1.2節(jié)中由力學性能試驗擬合的表達式給定。

圖5 算例模型

2.2 計算流程與方法

模擬金屬橡膠靜態(tài)受壓過程時，利用循環(huán)計算過程，將加載（或卸載）過程分為多個計算步，根據(jù)上一步的計算結果修正模型中的材料參數(shù)，更新彈性模量、泊松比和剪切模量，形成新的總剛度矩陣、質(zhì)量矩陣和載荷向量，并進行下一步循環(huán)計算，計算流程如圖6所示。采用APDL語言編寫計算程序，實現(xiàn)以上計算過程。

圖6 增量加載法計算流程

3 立方體件的受壓過程仿真

3.1 不同方向受壓過程

對金屬橡膠立方體件在成型方向和非成型方向上的受壓過程分別進行仿真計算，模擬其加載過程，并與試驗結果進行對比分析。

3.1.1 成型方向

模擬金屬橡膠立方體件在成型方向受載時，采用增量加載法對靜態(tài)加載過程進行模擬，加載至應變?yōu)?4%左右，計算輸出力-位移曲線，如圖7所示。從圖中可見，在24%應變范圍內(nèi)，仿真計算得到的力-位移曲線與試驗實測曲線吻合度高，計算結果較為精確。在金屬橡膠受壓過程中，仿真模擬的壓載荷-位移曲線與試驗實測曲線的誤差值呈現(xiàn)先增大后減小的趨勢。仿真模擬全程誤差在8%以內(nèi)。

圖7 成型方向計算結果與試驗結果對比

3.1.2 非成型方向

模擬金屬橡膠立方體件在非成型方向上加載時，采用增量加載法計算輸出的力-位移曲線如圖8所示。從圖中可見，在非成型方向加載至6%應變范圍內(nèi)，仿真計算結果與試驗所得曲線吻合度也較高。從誤差分布來看，在小載荷、小應變階段誤差較大，而當載荷增大后則誤差明顯減小。在仿真模擬的全過程中，誤差均在10%以內(nèi)。

圖8 非成型方向計算結果與試驗對比

3.2 約束狀態(tài)的影響

為研究約束狀態(tài)的影響，對立方體件在成型方向上的受壓過程進行深入研究，改變模型及試驗件的側向約束狀態(tài)，并對仿真和試驗結果進行對比分析。

3維立方體件在成型方向受壓時，承載面為各向同性平面，因此側面約束有3種狀態(tài)，即無側向約束、兩側約束（僅在1個非成型方向上約束）和四周約束（在2個非成型方向上均施加約束），如圖9所示。

圖9 金屬橡膠模型約束狀態(tài)（成型方向受壓）

對應3種約束狀態(tài)，在所建立模型上施加對應的位移約束，經(jīng)計算得到所施加的力和相應的位移曲線，如圖10所示。從圖中可見，在同樣大小的壓力下，無側向約束時得到的位移最大，兩側約束時略有減小，而當四周均被約束時則所得位移值明顯減小。在加載初期（位移小于3 mm時），3條曲線的差別不大，而隨著應變進一步增大，3條曲線的差距也逐漸增大。同時，當僅有1個非成型方向被約束時，力-位移曲線僅有極小幅度偏移，而當2個非成型方向均被約束后，曲線則有較為明顯的向上偏移，宏觀表現(xiàn)出剛度的非線性增大更加劇烈。

圖10 帶約束3D模型在成型方向的力-位移曲線

相同規(guī)格的試驗件在不同約束狀態(tài)下的試驗結果如圖11所示。從圖中可見，當金屬橡膠立方體件的側面存在約束時，其壓力-位移曲線發(fā)生偏移。約束狀態(tài)為兩側約束時，曲線相對無側向約束時偏移量并不大；而當立方體件的四周均被約束時，其壓力-位移曲線發(fā)生明顯偏移。此規(guī)律與仿真計算結果一致，驗證了仿真計算結果的準確性。

圖11 試驗實測帶約束金屬橡膠件的試驗力-位移曲線

為進一步定量評估仿真計算模型的精度，取約束狀態(tài)為四周約束的仿真與試驗曲線，對其剛度進行分析。二者的剛度對比及誤差分布如圖12所示。從圖中可見，在加載的初始階段，仿真模型的剛度低于試驗實測剛度，誤差值為25%左右。隨著加載應變的增加，二者的相對誤差略微增加，在應變值達到10%之后，仿真計算模型的相對誤差迅速減小。在20%應變范圍內(nèi)，采用有限元仿真計算模型得到的預估剛度值，其相對誤差在30%以下。

圖12 四周約束狀態(tài)下仿真與試驗所得剛度與誤差分布

4 環(huán)形件的徑向受壓過程仿真

圖13 金屬橡膠有限元模型加載

對金屬橡膠環(huán)形件進行有限元仿真計算時，在第2.1節(jié)所建立金屬橡膠環(huán)形件模型的基礎上，在其內(nèi)部建立剛性環(huán)，用以模擬實際工作狀態(tài)下的鼠籠軸頸或在試驗環(huán)境下的加載軸。如圖13所示，金屬橡膠環(huán)與剛性加載環(huán)之間建立接觸對單元，摩擦系數(shù)取0.15。在加載時，在金屬橡膠環(huán)外表面施加徑向約束，徑向壓載荷以集中力的形式施加在剛性環(huán)上，并由接觸對傳遞至金屬橡膠環(huán)。

按照金屬橡膠環(huán)形件在實際使用中的特點，將承載方向定義為非成型方向，軸向定義為成型方向，材料參數(shù)按橫向各向異性材料模型確定。金屬橡膠環(huán)在徑向加載時，經(jīng)有限元仿真計算和試驗測試獲得其壓力-位移曲線，如圖14所示。

圖14 金屬橡膠環(huán)徑向壓力-位移仿真與試驗結果對比

對圖14計算得到的數(shù)據(jù)進行處理，得到金屬橡膠環(huán)在8%應變范圍內(nèi)的剛度，如圖15所示。從圖中可見，當金屬橡膠環(huán)形件在徑向方向受壓時，其加載過程剛度隨應變的增大逐漸上升，但仿真計算所得剛度與試驗實測剛度相比偏小，相對誤差隨加載應變的增大而逐漸減小，在加載范圍內(nèi)均處于相對誤差25%之內(nèi)。仿真計算誤差的主要來源是，在計算中為提供求解效率和收斂性，并未考慮幾何變形非線性帶來的影響，因此有限元計算結果的變形值整體偏大，進而帶來剛度偏小。

5 結論

圖15 金屬橡膠環(huán)的徑向剛度及誤差分布

（1）在立方體件無約束自由狀態(tài)下，采用所建立的模型和計算方法得到的結果與試驗測試結果吻合。仿真模擬曲線與試驗實測曲線相比，無論是在成型方向上受壓還是非成型方向上受壓，相對誤差均在10%以內(nèi)。

（2）對于帶有側面約束的立方體件，約束會造成其加載過程中剛度增大，約束的自由度越多，則剛度的增加幅度越大。采用本文建立的模型和計算方法得到的結果與試驗測試結果相比，剛度的相對誤差在30%以內(nèi)。

（3）對于金屬橡膠環(huán)形件，在其徑向受壓時，采用本文建立的模型和計算方法得到的結果與試驗測試結果相比，剛度的相對誤差在25%以內(nèi)。

本文所建立的金屬橡膠有限元模型和仿真計算方法能充分反映金屬橡膠材料的橫向各向同性和應變相關性的非線性特征，適用于帶有復雜結構形式金屬橡膠的機械結構的靜力學計算，并具有較高的求解精度。在后續(xù)研究中，將進一步考慮金屬橡膠的阻尼特性，發(fā)展適用于動力學的仿真方法。

[1]Luis S A,Chirathadam T A.Metal mesh foil nearing:effect of motion amplitude,rotor speed,static load,and excitation frequency on force coefficients[J].Journal of Engineering for Gas Turbines and Power,2011,133(12):1201-1228.

[2]Ertas B H.Compliant hybrid journal bearings using integral wire mesh dampers[J].Journal of Engineering for Gas Turbines and Power,2009,131(2):1215-1226.

[3]Gadot B,Martinez O R,Roscoat S R D,et al.Entangled single-wire NiTi material:A porous metal with tunable superelastic and shape memory properties[J].Acta Materialia,2015,96:311-323.

[4]Li Y L,He Z B,Bai H B,et al.Advances in research and application of metal rubber[J].Ordnance Material Science and Engineering,2011,1:39.

[5]馬艷紅,陸宏偉,朱海雄,等.彈性環(huán)金屬橡膠支承結構剛度設計與試驗驗證[J].航空學報,2013,34(6):1301-1308.MA Yanhong,LU Hongwei,ZHU Haixiong,et al.Structural stiffness design and experimental evaluation of elastic ring metal rubber damper[J].Acta Aeronautica et Astronautica Sinica,2013,34(6):1301-1308.（in Chinese）

[6]洪杰,朱彬,馬艷紅,等.金屬橡膠顆粒靜態(tài)特性試驗[J].北京航空航天大學學報,2012,38(2):218-221.HONG Jie,ZHU Bin,MA Yanhong,et al.Experiments on static characteristics of metal rubber particles[J].Journal of Beijing University of Aeronautics and Astronautics,2012,38(2):218-221.（in Chinese）

[7]Liu P,He G,Wu L H.Uniaxial tensile stress-strain behavior of entangled steel wire material[J].Materials Science and Engineering A,2009,509:69-75.

[8]Zarzour M,Vance J.Experimental evaluation of a metal mesh bearing damper[J].Journal of Engineering for Gas Turbines and Power,2000,122(2):326-329.

[9]Zhang D Y,Scarpa F,Ma Y H,et al.Compression mechanics of nickel-based superalloy metal rubber[J].Materials Science and Engineering A,2013,580:305-312.

[10]Zhang D Y,Scarpa F,Ma Y H,et al.Dynamic mechanical behavior of nickel-based superalloy metal rubber[J].Materials and Design,2014,56:69-77.

[11]李宇燕,黃協(xié)清,宋凱.金屬橡膠非線性干摩擦副的接觸作用機理及其仿真結果分析[J].振動與沖擊,2011,30(7):77-81.LI Yuyan,HUANG Xieqing,SONG Kai.Contacting mechanism of nonlinear friction pair for metallic rubber and its simulation results[J].Journal of Vibration and Shock,2011,30(7):77-81.（in Chinese）

[12]于鋒禮.鋼絲網(wǎng)墊減振器的建模研究 [D].南京:南京航空航天大學,2008.YU Fengli.Study on modeling of steel-net pad damper[D].Nanjing:Nanjing University of Aeronautics and Astronautics,2008.（in Chinese）

[13]張玲凌,楊智春,孫浩.一類金屬橡膠阻尼器的建模與參數(shù)識別[J].機械科學與技術,2007,26(5):558-562.ZHANG Lingling,YANG Zhichun,SUN Hao.Modeling of a type of nonlinear metal rubber damper and its parameter identification[J].Mechanical Science and Technology for Aerospace Engineering,2007,26(5):558-562.（in Chinese）

[14]朱彬,馬艷紅,張大義,等.金屬橡膠遲滯特性本構模型研究[J].物理學報,2012,61(7):474-481.ZHU Bin,MA Yanhong,ZHANG Dayi,et al.A constitutive model of metal rubber based on hysteresis property[J].Acta Phys.Sin.，2012,61(7):474-481.（in Chinese）

[15]Ma Y H,Zhang Q C,Zhang D Y,et al.The mechanics of shape memory alloy metal rubber[J].Acta Materialia,2015,96:89-100.

[16]王亮.金屬橡膠密封系統(tǒng)模型建立及仿真分析 [D].黑龍江:哈爾濱工業(yè)大學,2010.WANG Liang.Metal rubber sealing system model establishment and simulation analysis[D].Heilongjiang:Harbin Institute of Technology,2010.（in Chinese）

[17]Ma Y H,Hong J.Study on metal rubber material’s characteristics of damping and sound absorption[R].ASME 2008-GT-50961.

[18]Tan Q B,He G.Stretching behaviors of entangled materials with spiral wire structure[J].Materials and Design,2013,46:61-65.

[19]Ma Y H,Gao D,Zhang D Y,et al.Compressive and dissipative behavior of metal rubber under constraints[J].Physica Status Solidi(b),2015,252(7):1675-1681.

Finite Element Simulation Method of Metal Rubber Based on Tested Parameter Identification

LUO Ma1,ZHANG Da-yi2,ZHU Hai-xiong2,MA Yan-hong2,3,HONG Jie2,3
（1.China Aviation Powerplant Research Institute,Zhuzhou Hunan 412002,China； 2.School of Energy and Power Engineering,Beihang University，Beijing 100191， China； 3.Collaborative Innovation Center of Advanced Aero-Engine,Beijing 100191， China）

During the metal rubber's(MR)engineering application,a numerical method with fast solving speed that is suitable for complex configurations is required.It is also hoped that the method could be combined with the commercial finite element software considering the engineers'requirements.Aiming at solving the finite element simulation problems of MR,the paper developed a parameter identification method for the transverse isotropy materials based on the tested data.The nonlinear mechanical characteristics of MR during compression process were calculated through the looping loading and model updating method.The results reveal that the relative error between predicted and tested result of cubic MR model without constraint is less than 10%and less than 30%for cubic MR under constraint.The constraint also leads to larger stiffness of MR.The relative error of ring MR model under radial loading is holded within 25%.

metal rubber;mechanical property;finite element;simulation;complex forced state

V 257

A

10.13477/j.cnki.aeroengine.2017.03.011

2016-05-19 基金項目：國家自然科學基金（51475021、51475023）資助

羅瑪（1978），男，高級工程師，主要從事航空發(fā)動機結構及性能設計方面工作；E-mail:rome158@163.com。

羅瑪，張大義，朱海雄，等.基于試驗數(shù)據(jù)識別的金屬橡膠有限元仿真方法[J].航空發(fā)動機，2017,43（3）：56-61.LUOMa，ZHANGDayi，ZHU Haixiong，et al.Finite element simulation method ofmetal rubber based on tested parameter identification[J].Aeroengine，2017,43（3）：56-61.

（編輯：栗樞）