本文主要依据个人之前学习和工作积累进行橡胶类超弹体材料本构模型在CAE仿真计算的技术简要整理，由于非此类科学技术计算专业工程人员，唯恐出错以致误导，诚恳大家辨别学习，但还是愿意以此种方式进行编写分享以能给那些曾经如我一样没人指导和参与培训学习的朋友。

 

很早之前就曾想过把CAE仿真中对于超弹体和粘弹性计算的准备和部分流程方法进行一个分类整理并分享出来。但自身仅做为机械工程普通从业人员，又不是专业的高分子材料仿真从事人员，底气始终是少了一点，即使曾经很长一段时间自学过相关理论和做过不少相关计算模型，所幸主导过百万费用级橡胶密封计算项目，但始终怕误导刚入门这方面仿真的朋友，迟迟不敢下笔。

几天前，仿真秀平台一句“如果说是真心的分享技术，那么读者自然会有自己的理解和判断，无须担心”，打消了我的顾虑。所以以这种心态利用闲散时间对超弹体和粘弹性计算的CAE仿真初步工作做一个整理，若是对部分读者有所帮助，那就是莫大的欣慰。本文分为两个篇幅，第一个篇幅进行橡胶类超弹体本构仿真计算的内容简述，另一个篇幅简述粘弹性仿真计算的准备工作。

 

超弹体仿真材料处理与本构拟合

1、橡胶超弹体材料处理

分析橡胶类仿真计算一般需要选择超弹体材料模型，超弹体材料模型假设材料是各向同性的、等温和弹性的，完全或接近不可压缩，是真实橡胶行为的理想化。

对于橡胶超弹体进行模拟仿真应该首要进行材料曲线的拟合工作，主流CAE仿真软件都提供了曲线拟合工具，可以帮助把实验数据转化成各种超弹模型能使用的应变能量密度函数系数。对于超弹体的试验数据种类可以选择图1中所示的多种或者至少一种，一般认为能够提供的数据种类越多，拟合的曲线越能表现真实橡胶特性，但对于以压缩为主的仿真计算项目，建议试验数据应该包括单轴压缩或等双轴拉伸。

图1

应该注意的是用于拟合曲线的测试的数据（除体积测试数据）需要工程应力-应变数据，体积测试数据需要真实应力-应变数据。这与金属非线性计算中塑性曲线拟合中收集数据方式不同，一般金属塑性曲线需要工程应力-应变转换为真实应力-应变。对于试验数据需要调整滞回和应力软化行为，采用稳定的曲线（应该偏移至零应力和零应变）用来进行曲线拟合，如图2所示。

图2

2、橡胶超弹体材料本构拟合

超弹体材料本构模型很多，如图3所示。

图3

同时同种本构模型又有不同的项数，例如常用的Mooney Rivlin 模型就有分为N=2，3，5和9项模型，可看作是多项式形式的特殊情形。

• 2项Mooney Rivlin模型相当于N=1的多项式形式，是最常用的模式之一。

• 3项Mooney Rivlin模型与N=2且 的多项式形式类似。

• 5项Mooney Rivlin模型相当于N=2多项式形式。

• 9项Mooney Rivlin模型相当于N=3的多项式形式。

所有的Mooney Rivlin模型，初始体积模量和初始剪切模量是

Mooney Rivlin 模型2项在拉伸应变大致90-100%时有效，但是2项MR不能很好描述压缩行为特性；必须保证提供足够的高阶项数据，捕捉任何拐点，5和9项可以用于应变高达100-200%的情况。但是当不能很好判断何种材料本构模型适配试验数据的状态下，仅个人习惯做法是对于不同材料本构模型与试验数据进行多次拟合尝试，选择试验数据贴合最高或最关心数据段计算状态吻合的材料本构模型，例如计算项目以压缩趋势为主，需要包括单轴压缩或者等双轴拉伸数据进行拟合，选择压缩趋势最吻合而不是剪切数据拟合吻合度高的材料本构，更好的建议应该了解不同材料本构模型计算局限适用应变范围，如图4部分总结。

图4

3、超弹体材料拟合流程

个人工作使用过的CAE软件中，Abaqus相对于ANSYS在超弹体材料本构拟合中可以一次性拟合出多个本构模型曲线以方便观察选择比较，限于个人原因依旧以自己较熟悉ANSYS平台对试验数据进行超弹体材料曲线拟合基本流程进行说明，具体流程如图5所示。

(1)【Engineering Data】工具箱中查找到Hyperelastic Experimental中Uniaxial Test Data双击或者右键点击Include Property，将选中项添加到材料属性窗口。可以对其他测试数据类型重复操作（例如双轴测试数据、剪切测试数据等）。

(2) 将测试数据粘贴到“Table of Properties”窗口。

(3) 从【Engineering Data】工具箱中，选择相应超弹性模型进行曲线拟合。

(4) 右键Curve Fitting并执行“Solve Curve Fit”，将运行最小二乘法曲线拟合程序，找到更好的系数组合。

(5) 最小二乘法拟合可以选择归一化或者绝对误差法。

(6) 右键Curve Fitting，执行命令Copy Calculated Values to Property，完成参数录入。

图5

 

计算模型与收敛性控制

1、计算几何选择

关于超弹体仿真计算，超弹体材料具有非常小的压缩性，假设完全不可压缩一般是非常好的近似，材料处理为几乎或完全不可压缩的选择取决于用户和可用的数据。超弹体的仿真计算是相对收敛困难的非线性材料计算类型，部分的计算项目收敛难度更是让工程师苦不堪言。

在进行处理计算几何选择方面，对于橡胶类超弹体计算，个人倾向于采用2D几何计算模型，尽力避免3D模型几何。

• 3D计算模型相对2D模型计算采用成数倍数量的单元节点，同时非线性计算的迭代求解和存储耗费，这都给计算资源和效率带来了挑战。

• 2D模型几何边线更方便控制尺寸份数以细化几何拐角，有效防止挤压或者拉伸时候的单元扭曲，提高收敛可能性。

• 更好处理橡胶大变形状态下，橡胶与金属件接触关系以及橡胶自接触；

• 规避开3D模型实体六面体单元的可能“刚硬”。

如图6和7所示是两类形状密封圈挤压状态模拟计算案例，模型计算都不可避免设置模型自接触、拐角单元细化以及单元形状控制以防止过度畸形。

图6

图7

2、不可压缩性的考虑事项

然而工程项目不能规避掉橡胶类计算都能很好简化成2D几何，仿真计算不可避免遇到3D橡胶超弹体计算模型。

这里仅代表个人讲下工作习惯，3D超弹体计算偏选择四面体单元而不用六面体单元，偏选择低阶单元而不用高阶单元，依靠增强应变单元解决剪切锁定（单元的位移场不能模拟由于弯曲而引起的剪切变形和弯曲变形，当单元长度与厚度的数量级相同或长度大于厚度时，此现象会更加严重），依靠采用U-P公式解决体积锁死（一般发生在完全积分单元中，体积锁死是指完全积分单元受到过度约束时的一种锁死现象，如果材料是不可压缩的或是近似不可压缩的，完全积分单元可能会变得特别刚硬而不会产生体积变形，即所谓“体积锁死“）。如果模型因为采用Mixed u-p而不能满足收敛，可以通过定义参数引入少量的不可压缩性。图8中计算的3D模型单元采用了四面体单元，而图9中管壁厚度考虑六面体单元计算并进行了单元选择控制。

图8

图9

3、非线性控制计算

橡胶超弹体材料计算属于难于收敛的材料非线性计算问题，应该进行比较细致周全的非线性计算控制，这是计算收敛的重要因素，但限于篇幅这里不能一一把非线性控制设置进行说明，读者可以参阅笔者在仿真秀之前编写的非线性控制文稿，此处不再累述，略。

4、自适应网格设置

即使考虑了上面的诸多计算因素，然而2D/3D超弹体仿真计算实际计算过程仍旧可能遇到因为单元畸形发散的问题。

这里简要提及一种自适应网格重构技术，可以作为非线性材料计算单元经历大变形、单元扭曲条件或提高求解精度的一种选择技术（但是不一定能成功），它基于一定准则能进行网格重画分，如图10所示的挤压计算，解决部分在大变形状态单元扭曲而无法满足收敛。

图10

主流的CAE软件一般都具备了非线性自适应网格的计算能力，自适应网格计算的搭建流程见图11所示，这里以个人较为熟悉的ANSYS来说一下部分了解的菜单，如图12所示。

图11

(1) 分析求解通用设置需求

需要定义非线性自适应网格可以实现的分析设置需求。

需要定义大变形控制打开Large Deflection=On，以及存储所有求解结果计算点Store Results At= All Time Points，部分读者也可以打开重启动设置以帮助发散后重新处理求解子步等。

图12

（2）定义自适应网格标准

自适应网格标准有3种准则，如图13所示。

图13

Energy：如果单元应变能大于等于组件中平均应变能乘以自定义参数，单元将会分割与重划。以获得高精度区域模拟结果。

Box：在指定区域内，进行分割和重划。

Mesh：有三种选项控制，其中

• Skewness：A larger value reduces the chance of a remesh to take place. The default value is 0.9.

• Jacobian Ratio：This value defines the maximum Jacobian ratio at element integration points. The value must be between -1 and 1. The default value is 0.1。

（理想值为1，大于0.7可接受，质量较好，小于0.5，准确性不能保证）

（3）Check At重划分点监测以及Time Range

重划分点监测能够选择Equally Spaced Points以及Specified Recurrence Rate来进行重化网格发生的次数与频率等，如图14所示。

时间范围控制包括Entire Load Step和Manual两种形式，当定义Manual时，能够定义起始时间和结束时间。另外，能够定义多组Nonlinear Adaptive Region，对不同时间内划分方法进行相应控制，如图14所示。

图14

（4）Nonlinear Adaptivity Remeshing Controls

当采用自适应网格技术的时候，在“Static Structural analyses→Analysis Settings”路径下将会出现Nonlinear Adaptivity Remeshing Controls自适应网格控制，能够对多个自适应网格项进行设置控制，如图15所示。限于篇幅和个人对这里技术理解水平较低，建议读者参阅相关帮助文件进行更深入的学习，例如其中Quality Improvement (NSL/GSR/RT)：控制重划网格去除扭曲，Refinement (NSL/GSR/RT)控制单元细化；例如Edge Angle（ESA）大于指定的edge angle临界值，将会进行分割，分割节点将会自动保留，并作为硬节点，等等。

图15

5、自适应网格应用案例

如图16所示是个人复现的ANSYS帮助中某超弹体计算的案例，案例采用自适应网格控制方法完成，可见即使压缩量如此极大，但是自适应网格能较好进行分析项目的完成，然后笔者提醒的是，自适应网格并不是压缩收敛的全部原因，这回到了文章前面部分介绍的M-R模型和Ogden模型应变计算的能力选择上，只有考虑超弹体本构能满足如此大的压缩应变的情况下，计算才能顺利收敛完成，请注意。

图16

 

流体渗压

 

图17

如图17所示，O型密封圈装入密封槽后，承受接触压缩应力后产生变形。介质压力不作用时，密封圈在自身的弹性作用下，对接触面产生一个预接触压力。充入压力介质后，密封圈发生位移，移向低压侧，弹性变形更大，填充和封闭密封间隙。作用在密封副耦合面的接触压力上升为

初装预接触应力能保证密封性能，具有自密封的作用。

这里最后介绍一种简易施加流体压力的方法，流体渗压技术作用的原因是因为结构受介质压力作用的过程中，接触介质的状态时刻改变。通过编写命令的方式进行起始点和闭合点的位置的搜索，一旦接触压力小于用户定义的穿透准则，接触面将会遭受流体的压力作用。而对于橡胶分析，需要进行多个载荷步和子步的求解，因此这种流体渗压的迭代过程在整个过程发生，技术细节内容如图18所示，案例复现效果如图19所示。

图18

图19