橡胶-帘线增强复合材料的本构模型及数值模拟

孙书蕾，陈文国

（1.贵州理工学院机械工程学院，贵州 贵阳 550003；2.曲靖师范学院信息工程学院，云南 曲靖 655011）

1 引言

橡胶-帘线复合材料在轮胎中有着广泛的应用。橡胶-帘线复合材料既有橡胶材料的特性，也有帘线的加强效果：橡胶是一种典型的超弹性材料，在帘线方向承担大的拉伸载荷。在轮胎的数值模拟中涉及到橡胶-帘线的分布以及材料体积的分配比例等，同时因为橡胶材料是典型的非线性材料以及接触非线性等多重非线性耦合分析，进行有限元的数值分析[1]相对比较困难，因此建立合理有效的橡胶-帘线等效本构模型尤为重要。研究表明，简单的按现有的复合材料层合板理论得到的结果与真实力学性能有很大差距。

文献[2]运用基于连续介质力学的方法，建立单位体积应变能本构模型，对一族纤维增强的橡胶帘线复合材料进行了数值模拟。文献[3]对橡胶-帘线复合材料单层板拉伸杨氏模量的预测方法进行了研究。文献[4]提出了一种建模方法。文献[5]对帘线断裂时的橡胶-帘线复合材料温度场数值模拟。文献[6]对超弹性基体中含有加强层的复合材料做了研究，通过均质化方法、在第三方向用解析的方法得到材料的等效力学性能。文献[7-8]把超弹性基体-纤维增强复合材料简化为横观各向同性的超弹性材料，运用数值模拟进行了分析，但是只是针对一族帘线加强的情况。

文献[9-12]运用基于连续介质力学的方法，提出一种新的非线性复合材料的本构模型，把人体椎间盘纤维环单位体积能量函数分为基质、纤维和基质两者的交互作用三部分，并通过分步实验确定本构模型的参数，实验和数值仿真的数据符合的很好。将这种方法应用到橡胶-帘线复合材料的数值模拟当中，通过引入帘线角度方向的结构张量从而把单位体积的应变能分解为便于识别的体积变形、等容变形和各向异性偏量变形三部分，从而建立了橡胶-帘线复合材料的一种新的超弹性本构模型，如图1 所示。并提出了识别其中材料参数的分步数值实验方案。进而通过数值算例，预测了帘线在不同排列方向条件下的力学性能。

图 1 轮胎帘线结构示意图Fig.1 Sketch Map of Tyre Cord Structure

2 橡胶-帘线本构模型

2.1 各向异性超弹性纤维增强本构模复合材料型

根据文献[13]的各向异性纤维增强复合材料的理论，单位体积的Helmholtz 应变能量函数依赖于变形梯度和纤维方向：

式中：C—右柯西-格林变形张量；A0=a0⊗a0—结构张量；a0—纤维初始方向的单位向量。

单位体积的Helmholtz 应变能量函数也可以写成关于C和a0的不变量的形式：

前三个不变量是各向同性的基质材料的贡献，后两个不变量与纤维材料的贡献相关。对于两族纤维的材料，如图1 所示。应变能的形式为：

式中：B0=b0⊗b0—关于纤维方向b0的结构张量。写成不变量的形式为：

变形梯度F可以分解为等容变形和体积变形，因此可以把应变能量函数分解为以下三个部分：

其中，α＝4，5 和 α0相关，α＝6，7 和b0相关，α=8 与两族纤维的方向相关。Wvol，Wiso，Wani分别代表体积变形，各向同性变形和各向异性变形。J是变形体积比为相应的修改的不变量：

为了便于参数拟合，根据文献[9]的分析，也可以把应变能分解为：

其中，WM想当于式中的相当于式中的

2.2 本构模型在轮胎中的应用

轮胎中有两族帘线，角度分别为30°、150°，如图1 所示。基质橡胶通常被看作是一种超弹性材料，在这里，采用比较典型的近似不可压缩的Mooney-Rivlin 模型：

纤维部分贡献的能量与两族帘线的角度相关，如图1 所示。

为方便起见，对于式中的Wani部分也即式中的WF+WFM部分，采用下面的简化形式：

本构方程确定以后，对应变能函数根据链导法则对右柯西-格林变量C求导数，就可以得到二阶Piola-Kirchhoff 应力张量：

对于近似不可压缩材料：

把带到可以得到：

对于所提出的单位体积的应变能量函数，可以推导出总的柯西应力表达式为：

式中：b=FFT—左柯西-格林张量—矩阵的迹；dev（·）=（·）-（1/3）［（·）：I］I，I—单位矩阵。

下面给出本构模型中参数的确定方案：

（1）拟合纯橡胶材料参数C10和C01。

（2）确定与体积模量有关的材料参数D。

（3）得到基质材料的参数以后，拟合沿帘线方向的单轴拉伸实验数据，得到C3和C4。

3 本构模型的参数拟合

3.1 橡胶部分的材料参数

运用三维数字散斑动态应变测量分析系统对纯橡胶部分进行单轴拉伸，如图2 所示。再运用ABAQUS 超弹性部分自带的数据拟合模块[14-15]，可以得到Mooney-Rivlin 的两个参数分别为：

C10=-0.048（MPa），C01=1.04（MPa）

图2 橡胶的单轴拉伸实验Fig.2 Uniaxial Tensile Test of Rubber

拟合曲线和实验结果对比，如图3 所示。

图3 橡胶单轴拉伸应力应变曲线图Fig.3 Stress-strain Curve of Rubber Under Uniaxial Tension

3.2 确定关于体积模量的材料参数D

设定橡胶的泊松比υ 为0.49，根据下式计算D[15]：

式中：μ0—初始剪切模量；

K—初始体积模量。

可以得到：

D=0.0198（MPa）-1

3.3 帘线部分的参数

得到橡胶材料部分的参数以后，通过拟合沿帘线方向的单轴拉伸实验数据[17]得到C3和C4，C3=C4=5997（MPa）。

拟合曲线与实验结果，如图4 所示。

图4 沿帘线方向单轴拉伸应力应变曲线图Fig.4 Curves of Uniaxial Tensile Stress-strain Curves along the Line of the Tire

4 数值算例分析

结合前面参数拟合得到的各项应变能的参数最终得到总的应变能表达式，再利用公式、就可以分别得到应力的表达式，在COSMOL 编制相应的自定义方程，模拟橡胶复合材料的单轴拉伸，如图5 所示。

图5 数值算例几何模型示意图Fig.5 Sketch Map of Geometric Model of Numerical Example

两族帘线对称布置，不同角度时的单轴拉伸得到的应力-应变曲线图，如图6 所示。

图6 帘线不同方向时的应力-应变曲线图Fig.6 Stress Strain Curve in Different Directions of Cord

从图6 中可以看出，不同的帘线角度时，对材料的力学性能影响很大。30°对称布置时得到的工程应力明显比60°对称布置时要高很多，并且随着应变的增大，两者之间的差距越来越大。而介于两者之间45°排列得到的工程应力则远远低于其他度角对称排列的到的结果，由此可见随着帘线角度的增大，材料强度随帘线角度的增加而减小，而在45°以后，则随帘线角度的增加而增加。这是由于橡胶材料的非线性本构关系，使得在不同帘线角度下，纤维和基体承担载荷的比例在发生变化，在帘线角度较小的情况下，载荷主要由帘线承受，帘线的承载能力明显高于基体。而在45°后，基体承载逐渐增强，所以应力值增加。对于这里的本构来说，45°是个转折点。这与文献[17]得到的结论是一致的。

5 结论

通过前面得到的本构方程和参数编制COMSOL 用户自定义方程，最后对不同帘线角度方向的轮胎材料进行单轴拉伸的数值模拟，结果表明不同的帘线角度时，对材料的力学性能影响很大，随着帘线角度的增大，材料强度随帘线角度的增加而减小，而在45°以后，则随帘线角度的增加而增加。

本构模型只用到了不变量I4和I6，并没有考虑到帘线和基质的相互作用，所以得到的结果不够精准，只能预测趋势。下一步的工作是将帘线和基质相互作用的的因素考虑进来，以更精确地预测橡胶-帘线材料的力学行为。