基于Exp-ln模型與廣義黏彈性理論的橡膠本構模型及其應用研究

仲健林，任杰,馬大為

(南京理工大學機械工程學院，南京210094)

摘要：提出了一種描述橡膠材料不同應變率下力學響應的黏超彈本構模型。首先，利用Instron實驗機和SHTB實驗裝置，開展橡膠材料單軸拉伸實驗；其次，結合Exp-ln超彈性本構模型和廣義黏彈性方法，建立了橡膠材料黏超彈本構模型；再次，推導本構模型三維增量格式，編寫了用戶子程序(VUMAT)，驗證了本構模型的一維和三維有效性；最后，建立橡膠底座沖擊附加載荷計算數(shù)值模型，并將沖擊附加載荷實驗與數(shù)值仿真進行對比驗證。結果表明：單軸拉伸實驗與數(shù)值解吻合較好，沖擊附加載荷實驗值與仿真值誤差約為7%，驗證了黏超彈本構模型的正確性。

關鍵詞：橡膠材料；熱黏超彈性；拉伸實驗；材料子程序(VUMAT);橡膠底座

中圖分類號:TJ762.1

文獻標志碼：A

DOI:10.13465/j.cnki.jvs.2015.19.024

Abstract:A visco-hyperelastic constitutive model was proposed to describe mechanical responses of rubber material under different strain rates. Firstly, the uniaxial tensile tests for rubber material were conducted with Instron testing machines and SHTB testing devices. Secondly, combined with Exp-ln hyperelastic constitutive model and the generalized viscoelastic method, the visco-hyperelastic constitutive model for rubber material was built. Thirdly, the three-dimensional incremental format of the visco-hyperelastic constitutive model was deduced, the user subroutine VUMAT was written, the one-dimensional effectiveness and three-dimensional one were verified. Finally, the impact additional load numerical simulation model of a rubber base was established, the results of impact additional load tests were compared with those of the numerical simulation. The results showed that the numerical results agree well with the uniaxial tensile test data, the error between the impact additional load test value and the simulation one is approximately 7%, the correctness of the visco-hyperelastic constitutive model is verified.

Constitutive model and its application for rubber material based on Exp-ln model and generalized viscoelastic theory

ZHONGJian-lin,RENJie,MADa-wei(School of Mechanical Engineering, Nanjing University of Science & Technology, Nanjing 210094, China)

Key words:rubber material; visco-hyperelastic constitutive model; tensile test; VUMAT; rubber base

橡膠是由重復單元(鏈節(jié))構成的長鏈分子，常用作橡膠底座等復合材料制品的基體材料[1]。從本質(zhì)上說橡膠是一種黏彈性材料，力學特性與應變率有關[2-3]。在準靜態(tài)條件下，橡膠的蠕變?nèi)崃炕蛩沙谀Ａ繛槌?shù)，橡膠表現(xiàn)為超彈特性，因此，橡膠材料的本構研究包括超彈性和黏彈性。Yang等[4-6]基于Yeoh函數(shù)提出了橡膠類材料的黏超彈本構模型并將其應用于沖擊仿真。Antoun[7-9]將本構方程表示為與坐標系無關的形式，提出了描述橡膠類材料黏彈性的本構模型。在小應變范圍內(nèi)(<30%)，上述本構模型的預測結果與實驗數(shù)據(jù)吻合較好，但由于本構建模過程中引入了由應變不變量組成的多項式，參數(shù)較多且物理含義不明確，可能會出現(xiàn)多解狀況[10]。因此，構建一種預測應變范圍大、參數(shù)少且物理含義明確的黏超彈本構模型其意義顯而易見。

文中利用Instron實驗機和SHTB實驗裝置，獲得了不同應變率下的橡膠材料拉伸應力-應變曲線；建立了橡膠材料黏超彈本構模型，超彈性本構采用參數(shù)較少且物理意義明確Exp-ln應變能函數(shù)描述，黏彈性本構建立在應變、時間效應解耦的基礎上，避免了引入應變不變量多項式；擬合了黏超彈本構相關參數(shù)，推導了本構模型的三維增量格式，編寫了本構模型的用戶子程序，將單軸拉伸、圓筒扭轉(zhuǎn)實驗分別與數(shù)值解進行對比驗證；應用該本構模型的用戶子程序完成橡膠底座沖擊附加載荷計算數(shù)值模型的建立，并將沖擊附加載荷實驗與數(shù)值仿真結果進行對比驗證。研究方法可為橡膠材料力學特性研究提供參考。

1不同應變率下的拉伸實驗

常見的獲取材料參數(shù)的實驗方法包括:單軸拉伸(壓縮)、雙軸拉伸和平面剪切等，單軸拉伸實驗方法簡單、滿足理論研究需要。針對橡膠材料，開展不同應變率下的單軸拉伸實驗，獲得應力-應變曲線，為黏超彈本構的參數(shù)擬合提供依據(jù)。

1.1橡膠材料單軸拉伸實驗

圖1　橡膠試件 Fig.1 rubbertest specimen

1.2實驗結果與分析

整理數(shù)據(jù)采集系統(tǒng)得到的信號，得到不同應變率下的工程應力P11-工程應變ξ11，見圖2。

圖2　不同應變率下拉伸實驗應力-應變曲線 Fig.2 Stress-strain curves of the tensile test under strain rate

由圖2可知，橡膠試件的應力隨應變的增大而增大；隨應變率的增大，橡膠試件的應力增大。橡膠材料的力學響應表現(xiàn)出明顯的應變率敏感性。

2橡膠材料黏超彈本構模型

橡膠材料黏超彈本構關系中包含率無關的超彈性應力σe和率相關的黏彈性應力σv。將超彈性本構模型與黏彈性本構模型進行組合，建立了橡膠材料黏超彈本構模型。

2.1Exp-ln超彈性本構模型

(1)

超彈性材料的應力應變關系可以用應變能函數(shù)W來度量，根據(jù)能量守恒可得，各向同性不可壓縮橡膠材料的本構方程為

(2)

常見的應變能函數(shù)W有：Yeoh應變能函數(shù)、Mooney-Rivilin應變能函數(shù)等。文獻[4]結合Mooney-Rivilin應變能函數(shù)給出橡膠材料的黏超彈本構模型，但是該模型需擬合的參數(shù)較多且預測的應變范圍有限。采用參數(shù)具有實際物理意義且適用應變范圍較廣的Exp-ln應變能函數(shù)，該應變能函數(shù)[11]為

(3)

(4)

(5)

(6)

2.2基于廣義黏彈性方法的黏超彈本構模型

目前，常見的黏彈性本構描述是將其表示為與坐標系無關的形式，使用應變不變量的多項式表示應變歷史對應力影響的張量泛函，多項式參數(shù)較多，物理意義不明確。本文結合Gamonpilas[12]給出的廣義黏彈性本構模型的描述方法，將應變ξ和時間t對應力的影響效應進行解耦，避免了引入多項式，黏超彈本構模型為

σ(ξ,t)=σ0(ξ)g(t)=σe+σv=

(7)

積分形式的本構方程能直接反映黏彈性材料的記憶特性，結合Leaderman卷積積分[13]和式(7)，橡膠材料黏超彈性本構為

(8)

將式(8)分解為超彈性應力部分和黏彈性應力部分，三維全量格式黏超彈本構模型為

σ=σe+σv=σ0(ξ)g∞(ξ)+

(9)

(10)

假設加載前材料的應力狀態(tài)不影響加載后的應力狀態(tài)，即時間積分下限為零，為減少模型參數(shù)數(shù)量，取一個松弛時間(N=1)，聯(lián)立式(9)、(10)可得

(11)

(12)

(13)

(14)

b′ln(λ2+2λ-1-2)]

(15)

通過式(15)結合動態(tài)單軸拉伸實驗數(shù)據(jù)，在求得兩個超彈性參數(shù)a,b的值的基礎上，利用最小二乘法擬合確定三個黏彈性參數(shù)A′,a′,b′和松弛時間θ1的值。文中橡膠的應變率跨越3個數(shù)量級，而橡膠的力學行為與高應變率、低應變率有關。本研究中取N=1，僅通過一個松弛時間來描述橡膠材料在整個應變率范圍的力學行為，可能會造成描述的誤差，有待進一步考證。

3本構模型的數(shù)值驗證

3.1黏超彈本構的參數(shù)擬合

將橡膠材料單軸拉伸實驗數(shù)據(jù)代入黏超彈本構模型中，進行最小二乘法擬合，即可獲得本構模型的6個待定參數(shù)，參數(shù)擬合步驟為

步驟1 根據(jù)準靜態(tài)加載條件下橡膠材料單軸拉伸實驗數(shù)據(jù)，結合式(6)并采用最小二乘法確定兩個超彈性參數(shù)a,b的值。

步驟2 將a,b代入式(15)中，選取應變率為1000s-1的曲線進行非線性擬合，確定三個黏彈性參數(shù)A′,a′,b′和松弛時間θ1的值。至此橡膠材料黏超彈本構模型中所有參數(shù)全部求得，見表1。

表1　黏超彈本構模型參數(shù)值

3.2本構模型三維增量格式與VUMAT開發(fā)

為進行VUMAT子程序開發(fā)，需對式(9)三維全量格式本構進行有限元離散，建立三維應力增量與應變增量之間的關系。

(17)

(18)

又

(19)

由式(17)得

dS=DijkldE

(21)

結合式(18)和式(20)，即可給出

Dijkl=4[We11δjiδlk+2We12δjigkl+3We13δjigkmgml+

2We21gijδlk+4We22gijgkl+6We23gijgkmgml+3We31gimgmjδlk+

6We32gimgmjδlk+6We32gimgmjgkl+9We33gimgmjgkngnl+

2We2δkiδlj+3We3(δkiglj+gikδli)]

(22)

將式(3)應變能函數(shù)We代入式(22)求得切線本構張量Dijklei，通過編制VUMAT子程序獲得Kirchhoff應力張量S和Green應變張量E的更新，對于不可壓縮橡膠有σe=FSFT、E=1/2(FTF-ξ)，即可獲得超彈本構模型σe的應力更新。

(23)

(24)

對于式(24)采用分步積分，第一部分積分起止時間為(0，t)，第二部分積分起止時間為(t，Δt)，當Δt很小時，近似認為t和Δt時刻的應變率相等，則式(24)可化簡為

(25)

至此，將超彈本構模型、黏彈本構模型的應力更新進行疊加，即獲得了進行子程序開發(fā)所需的本構方程三維增量格式，子程序開發(fā)的任務就是在材料積分點上完成上式的應力計算，并對主程序進行應力更新。針對式(21)、式(25)和VUMAT所提供的數(shù)據(jù)接口，制定黏超彈本構模型子程序開發(fā)實施步驟為：

步驟1從ABAQUS主程序中讀入增量步開始時刻材料屬性：abA′a′b′θ1；變形梯度Ft、Ft+Δt；時間增量Δt等；

步驟2從用戶自定義狀態(tài)變量載入Green應變張量E,根據(jù)E=1/2(FTF-ξ)求解Green應變增量dE；

步驟3由式(21)計算Kirchhoff應力增量dS，根據(jù)σe=FSFT求解獲得超彈本構模型σe的應力更新；

步驟4由式(25)計算獲得粘彈本構模型σv應力更新

步驟5將超彈本構模型σe的應力更新和粘彈本構模型σv應力更新進行疊加，并返回增量步結束時刻的應力；

步驟6更新狀態(tài)變量。

根據(jù)以上步驟，用FORTRAN語言編寫子程序，從而將橡膠黏超彈本構模型嵌入到ABAQUS中。

3.3本構模型的一維有效性驗證

為驗證本構模型的一維有效性，將3.2節(jié)中開發(fā)的橡膠黏超彈本構模型的VUMAT子程序用于模擬單軸拉伸實驗。見圖3橡膠試件有限元模型所示，采用C3D8R單元進行網(wǎng)格劃分，共計1866個單元，沖擊加載以實驗實測的入射桿應力波加載與入射桿端部，材料性能參數(shù)取表1中的數(shù)據(jù)。

圖3　橡膠試件有限元模型 Fig.3 Finite element model of the rubber test specimen

通過數(shù)值模擬，獲得了不同應變率下單軸拉伸實驗的數(shù)值解。將實驗數(shù)據(jù)與黏超彈本構的數(shù)值解進行對比驗證。給出不同應變率下的工程應力P11-工程應變ξ11對比曲線。見圖4，在工程應變ξ11<0.9范圍內(nèi)，工程應力數(shù)值解與單軸拉伸實驗解符合較好，表明文中提出的黏超彈本構模型能夠在較大范圍內(nèi)預測橡膠材料的單軸加載方向力學響應，驗證了黏超彈本構模型的一維有效性。

圖4　不同應變率下應力-應變曲線對比 Fig.4 Stress-strain curves under different strain rate

3.4本構模型的三維有效性驗證

為驗證本構模型的三維有效性，將3.2節(jié)中開發(fā)的橡膠黏超彈本構模型的VUMAT子程序用于模擬文獻[16]中橡膠襯套內(nèi)表面扭轉(zhuǎn)過程。見圖6(a)橡膠襯套有限元模型所示，采用C3D8R單元進行網(wǎng)格劃分，共計5760個單元，橡膠襯套內(nèi)徑9.85mm，外徑18.2mm，軸向長度60mm。橡膠襯套外圈固定，分別以30°/s、15°/s、7.5°/s、3.75°/s、1.875°/s的加載率在內(nèi)圈上施加15°的角位移，然后保持角位移不變。橡膠黏超彈本構模型的超彈性參數(shù)利用文獻[16]中的單軸靜態(tài)拉伸實驗數(shù)據(jù)進行擬合，黏彈性參數(shù)和松弛時間則利用文獻[16]中加載率為0.25mm/s的位移加載實驗數(shù)據(jù)進行擬合，獲得文獻[16]中橡膠材料黏超彈本構模型的參數(shù)擬合值，見表2。

表2　文獻[16]中黏超彈本構模型參數(shù)值

圖6　橡膠襯套有限元模型和松弛狀態(tài)應力云圖 Fig.6 Finite element model bushing and stress nephogram in relaxed state of the rubber

圖7　橡膠襯套內(nèi)表面扭矩時程曲線數(shù)值解與實驗解 [16] Fig.7 Numerical and experimental [16] results of the moment history curve for the inner surface of the rubber bushing

通過數(shù)值計算，獲得橡膠襯套松弛狀態(tài)下的應力云圖和內(nèi)表面的扭矩時程響應。不同角位移加載率下橡膠襯套松弛狀態(tài)下的應力云圖均見圖6(b)，這是由于橡膠材料的記憶衰退特性，在角位移加載保持恒定一段時間后，不同加載率下的橡膠襯套進入同一穩(wěn)定狀態(tài)；見圖7，將橡膠襯套內(nèi)表面受到的扭矩時程響應數(shù)值解和文獻[16]中實驗值進行對比，扭矩時程響應數(shù)值解和實驗解符合較好，驗證了橡膠黏超彈本構模型的三維有效性。

4本構模型的應用

為表明黏超彈本構模型描述高應變率效應的實用性，建立兩種選用不同本構模型(M-R超彈性本構、黏超彈性本構)的橡膠底座沖擊附加載荷計算模型，并將沖擊附加載荷實驗值與數(shù)值仿真進行對比驗證。

4.1沖擊附加載荷計算數(shù)值模型

橡膠底座常用于文獻[1]中所述的導彈懸垂彈射系統(tǒng)。見圖8中橡膠底座沖擊附加載荷計算數(shù)值模型所示，彈射沖擊瞬間，高壓氣體充入橡膠底座內(nèi)部，底座發(fā)生膨脹變形，最后與地面接觸，釋放壓力達到平衡狀態(tài)，作用在底座上的氣體壓力不能全部傳遞到地面，底座對初容室底部形成“附加載荷”。

圖8　沖擊附加載荷計算數(shù)值模型 Fig.8 Numerical calculation model for the impact additional load

橡膠底座由橡膠層和簾子布多層粘接熱壓而成，采用S4R單元對底座結構進行離散。采用Rebar單元模擬簾子布增強相；為對比應變率效應，分別采用M-R超彈性本構模型和本文提出的黏超彈本構模型用戶子程序來模擬橡膠材料；初容室和地面分別賦予相應的材料屬性。在數(shù)值模型中定義橡膠底座上端面和初容室底部固連，底座底面與地面之間建立接觸關系，在初容室與底座的連接界面上提取沖擊附加載荷作用力。

4.2附加載荷實驗與分析

見圖9，附加載荷測量實驗平臺由底座、初容室、測力傳感器、支撐臺架等部件構成。通過進氣閥在底座內(nèi)部進行加壓，底座膨脹與地面接觸后對初容室產(chǎn)生附加載荷。沖擊附加載荷通過初容室上端面?zhèn)鬟f至支撐臺架，利用周向均布的三個測力傳感器測量初容室上端面與支撐臺架之間的作用力，沖擊附加載荷值遠大于初容室的自身重力，該作用力可視為沖擊附加載荷。

圖9　實驗平臺示意圖 Fig.9 Schematic diagram of experimental platform

針對實驗模型，將三個測點的作用力曲線疊加得到?jīng)_擊附加載荷曲線，在橡膠材料參數(shù)不同的兩類中數(shù)值模型中分別提取附加載荷。對計算結果進行無量綱化處理，將附加載荷除以底座內(nèi)部壓強峰值p與底座上端面的面積S的乘積，圖10給出了沖擊附加載荷隨無量綱化時間變化曲線。

圖10　沖擊附加載荷對比 Fig.10 Comparison of additional impact load

由圖7可知，兩類橡膠材料參數(shù)不同的數(shù)值模型附加載荷數(shù)值解與實驗值變化規(guī)律一致，但M-R超彈性本構模型的計算結果與實驗存在明顯的相位差，而文中提出的黏超彈本構模型計算結果與實驗符合較好，附加載荷峰值的相對誤差約為7%，滿足工程精度要求，表明了黏超彈本構模型描述高應變率效應的實用性。

5結論

本文結合Exp-ln超彈性本構模型、廣義黏彈性方法，建立了橡膠材料黏超彈本構模型；通過用戶子程序驗證了其有效性；將黏超彈本構模型應用于橡膠底座的建模。得到了一些結論：

(1) 隨應變率的增大，橡膠材料單軸拉伸的應力增大，橡膠材料表現(xiàn)出應變率硬化效應。橡膠材料的力學響應與應變率密切相關。

(2) 黏超彈本構模型的超彈性部分采用參數(shù)少且物理意義明確的Exp-ln本構模型，黏彈性部分則建立在應變、時間效應解耦的基礎上。本文提出的黏超彈本構模型預測應變范圍大、參數(shù)少且物理意義明確。

(3) 在工程應變ξ11<0.9范圍內(nèi)，工程應力數(shù)值解與單軸拉伸實驗解符合較好；橡膠襯套松弛狀態(tài)下的應力云圖和內(nèi)表面的扭矩時程響應與文獻中實驗符合較好，驗證了文中提出的黏超彈本構模型及其用戶子程序的一維和三維有效性。

(4)相較于傳統(tǒng)的超彈性本構模型，采用黏超彈性本構模型的底座數(shù)值模型仿真結果與沖擊附加載荷的實驗曲線符合較好，峰值誤差約為7%，表明了文中黏超彈本構模型及用戶子程序描述高應變率效應的實用性。

文中研究方法和結論能夠為橡膠材料的力學響應預測提供技術支撐。

參考文獻

[1]仲健林,任杰,馬大為,等. 基于細觀力學精確建模方法的自適應底座力學性能研究[J]. 固體火箭技術,2014,3:400-407.

ZHONG Jian-lin,REN Jie,MA Da-wei,et al. Mechanical properties research on adaptive base based on micromechanics accurate modeling method[J]. Journal of Solid Rocket Technology,2014,3: 400-407.

[2]Drozdov A D, deC Christiansen J. Consti-tutive equations in finite elasticity of swollen elastomers [J]. International Journal of Solids and Structures, 2013, 50(9):1494-1504.

[3]Kovalev V A, Radaev Y N. Three-dimensional constitutive relations of ideal plasticity and the flow on the coulomb-trescaprism edge[J]. Mechanics of Solids, 2010, 45 (2):295-308.

[4]Yang L M, Shim V P W, Lim C T. A visco-hyperelastic approach to modeling the constitutive behavior of rubber[J]. International Journal of Impact Engineering, 2000, 24(6): 545-560.

[5]Yang L M, Shim V P W. A visco-hyperelastic constitutive description of elastomeric foam[J]. International Journal of Impact Engineering, 2004, 30(8):1099-1110.

[6]周相榮,王強,王寶珍,等.一種基于Yeoh函數(shù)的非線性黏超彈本構模型及其在沖擊仿真中的應用[J].振動與沖擊,2007,26(5):33-37.

ZHOU Xinag-rong, WANG Qiang, WANG Bao-zhen, et al. A nonlinear visco-hyperelastic constitutive model based on Yeoh strain energy function with its application to impact simulation[J]. Journal of Vibration and Shock, 2007,26(5):33-37.

[7]Pouriayevali H, Guo Y B, Shim V P W. A viscohyperelastic constitutive description of elastomer behavior at high strain rates [J]. Procedia Engineering, 2011, 10: 2280-2285.

[8]常新龍,賴建偉,張曉軍,等.HTPB推進劑高應變率黏彈性本構模型研究[J].推進技術, 2014, 35(1):123-127.

CHANG Xin-long, LAI Jian-wei, ZHANG Xiao-jun, et al.High strain rate visco-elastic constitutive model for HTPB propellant[J]. Journal of Propulsion Technology, 2014,35(1):123-127.

[9]Chen X H. Nonlinear electro-thermo-visco-elasticity [J]. Acta Mechanica,2010,211 (1/2): 49-59.

[10]胡少青，鞠玉濤，常武軍，等. NEPE固體推進劑粘一超彈性本構模型研究[J].兵工學報,2013,34(2):168-173.

HU Shao-qin, JU Yu-tao,CHANG Wu-jun, et al. A visco-hyperelastic constitutive behaviour of NEPE propellant[J]. Acta ArmamentarII, 2013,34 (2):168-173.

[11]Khajehsaeid H, Arghavani J, Naghdabadi R. A hyperelastic constitutive model for rubber-like materials[J]. European Journal of Mech-anics/A Solids, 2013,38: 144-151.

[12]Gamonpilas C,McCuiston R. A nonlinear viscoelastic material constitutive model for polyurea[J]. Polymer,2012,53(17):3655-3658.

[13]Li Shao-hua, Yang shao-pu,Chen Li-qun, et al.A nonlinear vehicle-road coupled model for dynamics research[J].Journal of Computational and Nonlinear Dynamics,2013, 8(2): 1-14.

[14]Anani Y, Alizadeh Y. Visco-hyperelastic consti-tutive law for modeling of foam’s behavior[J]. Materials and Design, 2011,32:2940-2948.

[15]張偉, 魏剛, 肖新科,等.2A12鋁合金本構關系和失效模型[J].兵工學報,2013,34(3):276-282.

ZHANG Wei, WEI Gang, XIAO Xin-ke, et al. Constitutive relation and fracture criterion of 2A12 aluminum alloy[J]. Acta ArmamentarII, 2013,34(3):276-282.

[16]Kadlowec J, Wineman A, Hulbert G. Elastomer bushing response: experiments and finite element modeling[J]. Acta Mechanica, 2003, 163(1/2): 25-38.