本发明涉及材料分析技术领域，特别涉及一种橡胶材料应力松弛及损伤效应的数值预测方法。
背景技术：
建立准确合理的本构模型是进行橡胶应力松弛-损伤特性分析的基础，所建立的本构模型将决定橡胶数值模拟能否全面表征材料的松弛行为及不可恢复损伤变形效应。此外，将所建立本构模型与现有数值模拟技术结合，构建准确的数值分析方法对橡胶制品的松弛性能预测尤为关键。
技术实现要素：
为解决上述技术问题，本发明提供了一种橡胶材料应力松弛及损伤效应的数值预测方法，以达到可对不同加载水平的橡胶材料应力松弛及损伤残余特性进行全面分析的目的。
为达到上述目的，本发明的技术方案如下：
一种橡胶材料应力松弛及损伤效应的数值预测方法，包括如下步骤：
第一步：橡胶材料时变松弛-损伤本构模型的建立；
第二步：橡胶材料的力学试验测试；
第三步：将时变松弛-损伤本构模型应用于有限元软件，构建时变松弛-损伤本构关系与有限元数值模拟技术的相容性仿真计算方法；
第四步：根据第二步获得的橡胶材料的力学试验数据，利用第三步建立的相容性仿真计算方法对橡胶材料时变松弛-损伤本构模型的关键参数进行拟合识别；
第五步：根据建立的时变松弛-损伤本构模型及识别的参数，通过第三步建立的相容性仿真计算方法，对橡胶材料的时变松弛行为及损伤变形数值进行模拟预测。
上述方案中，所述第一步的具体方法如下：
其中，
f是屈服函数，k为强化函数，为等效塑性应变，为等效塑性应变的增量，sij为应力偏张量，为等效屈服应力；橡胶材料塑性流动过程中的损伤可以分阶段表示，为不同载荷阶段的等效塑性应变，为不同载荷阶段的等效应力，当n＝0时，可得初始等效屈服应力
上述方案中，所述第二步的具体方法如下：
(1)开展橡胶试样在不同变形模式下的超弹力学试验，分别进行单轴拉伸力学测试、平面拉伸力学测试、等双轴拉伸力学测试及单轴压缩力学测试；
(2)采用与超弹力学试验测试相同批次的橡胶材料制成应力松弛试样，开展橡胶材料在不同加载位移水平的长时应力松弛试验；
(3)对橡胶试样在应力松弛卸载后的损伤残余变形进行测量，获取不同加载位移水平下的残余变形数据。
上述方案中，所述第三步的具体方法如下：
(1)在有限元软件abaqus中，建立橡胶试样的物理模型；
(2)对abaqus进行二次开发，基于fortran语言、uhyper子程序及abaqus脚本文件，将提出的时变松弛-损伤本构模型嵌入至商业有限元软件中，实现材料模型的编写入库；
(3)按照试验工况，对橡胶试样模型进行材料模型赋予、网格划分、边界条件和载荷条件施加、接触关系及求解计算设置，实现橡胶试样的有限元建模与分析。
上述方案中，所述第四步的具体方法如下：
(1)根据已获得的橡胶材料在单轴拉伸、平面拉伸、等双轴拉伸及单轴压缩变形模式下的超弹试验数据，采用最小二乘法对时变松弛-损伤本构模型中的超弹材料参数c01、c10和d10进行拟合识别；
(2)根据已获得的橡胶材料在不同加载位移水平的应力松弛试验数据，对时变松弛-损伤本构模型的损伤塑性应变材料参数进行拟合识别；
(3)在已识别超弹材料参数和部分损伤塑性应变材料参数的基础上，根据已获得的橡胶材料在不同加载位移水平的应力松弛试验数据，对时变松弛-损伤本构模型的松弛材料参数及损伤塑性材料参数进行拟合识别；具体包括：利用软件isight集成第三步中建立的橡胶材料相容性仿真计算方法，以时变松弛-损伤本构模型中的松弛材料参数k和r，和损伤塑性应力参数作为优化设计变量，以松弛试验结果与仿真结果的差值作为优化目标函数，采用多岛遗传优化算法开展材料参数拟合识别。
通过上述技术方案，本发明提供的一种橡胶材料应力松弛及损伤效应的数值预测方法具有如下有益效果：
(2)本发明所提供的橡胶材料应力松弛数值模拟方法具有建模与分析过程明确、预测精度高等优点，可对不同加载水平的橡胶材料应力松弛及损伤残余特性进行全面分析.
(3)本发明所提出的时变松弛-损伤本构模型的参数识别方法，具有易理解、易操作、拟合精度高的优点。
(4)本发明为橡胶材料及制品的实际工程力学性能分析提供了良好的预测模拟方法。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍。
图1为本发明橡胶材料应力松弛及损伤效应的数值预测方法流程图；
图2a-2d为丁腈橡胶不同变形模式下的试验应力应变曲线图，其中，图2a为单轴拉伸试验的应力应变曲线图，图2b为平面拉伸试验的应力应变曲线图，图2c为等双轴拉伸试验应力应变曲线图，图2d为单轴压缩试验的应力应变曲线图；
图3a-3d为氢化丁腈橡胶不同变形模式下的试验应力应变曲线图，其中，图3a为单轴拉伸试验的应力应变曲线图，图3b为平面拉伸试验的应力应变曲线图，图3c为等双轴拉伸试验应力应变曲线图，图3d为单轴压缩试验的应力应变曲线图；
图5为有限元模型图；
图6a-6h为不同拉伸试验和单轴压缩试验数据与数值拟合的结果对比图，其中，图6a，图6b，图6c，图6d分别为丁腈橡胶单轴拉伸、平面拉伸、等双轴拉伸、单轴压缩过程中的应力应变曲线对比图；图6e，图6f，图6g，图6h分别为氢化丁腈橡胶单轴拉伸、平面拉伸、等双轴拉伸、单轴压缩过程中的应力-应变曲线对比图；
图7为多岛遗传算法(miga)识别本构模型参数的流程图；
图8a-8b为应力松弛过程中应力衰减的数值模拟与试验测试结果对比图，其中，图8a为丁腈橡胶时变应力衰减曲线对比图，图8b为氢化丁腈橡胶时变应力衰减曲线对比图；
图9a-9b为应力松弛结束后损伤残余变形的数值模拟与试验测试结果对比图，其中，图9a为丁腈橡胶的损伤残余变形数值对比图，图9b为氢化丁腈橡胶的损伤残余变形数值对比图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述。
本发明提供了一种橡胶材料应力松弛及损伤效应的数值预测方法，如图1所示，具体实施例如下：
为验证本发明所构建的本构模型及数值模拟方法的可靠性，本实施例对两种不同类型的橡胶材料(丁腈橡胶和氢化丁腈橡胶)进行松弛特性数值预测。具体实施步骤如下：
第一步：时变松弛-损伤本构模型的提出
在工程应用中，橡胶往往处于大变形状态，其力学特征可采用基于应变能函数的超弹本构模型来模拟，基于应变不变量表示的连续介质唯象模型可写成，
式中，为材料加载响应的能量密度偏量函数，w(j)应变能密度的体积能函数。
可进一步分解成即第一和第二偏应变不变量,它们的定义是，
式中，是偏延伸变量。j是橡胶变形后与变形前的弹性体积比，即：
式中，b为左柯西格林变形张量，f为变形梯度张量。
根据上述公式可推导出不同的超弹本构模型，其中，基于多项式形式的超弹本构模型被广泛应用于工业界，可预测橡胶在外荷载激励下的力学响应。为了方便工程应用，本实施例中选定n＝2的多项式形式超弹本构模型(即目前广泛应用的mooney-rivlin模型)，其表达式为，
式中，σij为橡胶材料的应变能应力分量，为超弹与松弛应变增量。
材料的损伤塑性准则反映了材料产生塑性流动时的应力水平，可以用屈服函数来表示，如下公式所示：
式中，f是屈服函数，k为强化函数，为等效塑性应变，为初始等效屈服应力。
当f＜k时，橡胶材料处于时变超弹松弛状态，此时只产生超弹松弛变形，此时，
式中，dεij为总应变增量。
当f＝k时，橡胶材料进入超弹松弛-损伤塑性状态，即初始屈服后，橡胶材料的变形由超弹松弛变形与损伤塑性变形共同组成，具体如下，
式中，为损伤塑性应变增量。
进入超弹松弛-损伤塑性阶段后，随着加载位移的增加，损伤塑性变形也将进一步增加，而确定流动规则是确定残余变形变化趋势的前提条件。流动规则的增量形式如下，
式中，为等效塑性应变的增量，sij为应力偏张量，为等效应力。橡胶材料塑性流动过程中的损伤可以分阶段表示，为不同载荷阶段的等效塑性应变，为不同载荷阶段的等效应力，当n＝0时，可得初始等效屈服应力
第二步：静态拉伸试验、应力松弛试验
(1)测试丁腈橡胶及氢化丁腈橡胶在单轴拉伸、等双轴拉伸、平面拉伸及单轴压缩变形下的应力-应变关系，测试方案为：首先将橡胶材料裁成标准哑铃型试样，进行准静态单轴拉伸力学测试；其次，对平面状橡胶材料进行平面拉伸力学测试；再次，对橡胶材料进行等双轴拉伸力学测试；最后，对橡胶圆柱试样进行单轴压缩力学测试。
针对不同变形模式下的超弹力学试验，分别将橡胶试样安装在相应的试验机上，调整好试验机参数，进行测试。丁腈橡胶在不同变形模式下的应力-应变曲线如图2a-图2d所示，氢化丁腈橡胶在不同变形模式下的应力-应变曲线如图3a-图3d所示。
(2)丁腈橡胶及氢化丁腈橡胶的应力松弛试验
第三步：时变松弛-损伤本构关系与有限元数值模拟的相容性计算方法构建
通过二次开发及编译，将本发明所提供的时变松弛-损伤本构模型应用于有限元商业软件，构建时变松弛本构关系与数值模拟技术的相容性仿真计算方法。
(1)在商业有限元软件abaqus中建立橡胶试样的物理模型。参考试验样件为圆柱形试样(直径d为10.0mm，高度h为10.0mm)，建立橡胶试样的轴对称几何模型。
(2)对abaqus进行二次开发，基于fortran语言及uhyper子程序，将提出的时变松弛-损伤本构模型嵌入商业有限元分析软件中，实现材料模型的编写入库。将构建的应变能函数写入uhyper子程序，定义超弹本构模型中的第一应变不变量、第二应变不变量、弹性体积比、不可压缩性、超弹参数、松弛参数以及各变量之间的数值关系。各变量之间的计算关系如下，
式中，u为uhyper子程序定义的应变能函数，其中，松弛损伤部分作适当简化，得到如下关系，
损伤本构模型可通过abaqus脚本文件进行输入。
(3)按照试验工况，对橡胶试样模型进行材料模型赋予、网格划分、边界条件和载荷条件施加、接触关系及求解计算设置等，实现橡胶试样的有限元建模与分析。
对于本发明所提出的时变松弛-损伤本构模型，通过abaqus脚本文件将塑性参数等材料属性赋值给有限元模型，通过uhyper自定义程序将超弹和松弛等材料参数赋值给有限元模型。根据模型的对称性，采用四节点轴对称时变耦合四边形、双线性位移杂交单元来划分有限元模型的网格。橡胶试样模型的上、下表面分别与上、下解析刚体接触，位移载荷施加于解析刚体。构建的橡胶试样有限元模型如图5所示。
第四步：橡胶本构模型的关键参数识别
(1)辨识时变应力松弛-损伤本构模型的超弹材料参数
将单轴拉伸、平面拉伸、等双轴拉伸的试验数据导入abaqus，采用最小二乘法对时变松弛-损伤本构模型中的超弹材料参数(c01、c10和d10)进行数值拟合。对于每组应力应变试验数据，均使得相对误差ee最小，
式中，ee为相对误差，titest为试验应力结果，tith为按照mooney-rivilin超弹本构关系在与试验数据相同应变下对应的应力仿真结果。
将单轴压缩变形的超弹试验数据与上述数值拟合的结果进行对比，确定橡胶材料的超弹参数。丁腈材料参数：c10＝1414930pa，c01＝-169289pa，d10＝1.0×10-8；氢化丁腈材料参数：c10＝2030774pa，c01＝-870827pa，d10＝1.0×10-8。关于不同变形模式的超弹力学特性，丁腈试验与数值拟合的对比结果如图6a-6d所示，氢化丁腈试验与数值拟合的对比结果如图6e-6h所示。
(2)辨识时变应力松弛-损伤本构模型的损伤材料参数
第二步的试验中已记录橡胶试样松弛前后的尺寸数据，据此可确定时变松弛-损伤本构模型的损伤塑性应变材料参数丁腈橡胶的试验数据见表1，氢化丁腈橡胶的试验数据见表2。
表1丁腈材料松弛前后的厚度变形结果
表2氢化丁腈橡胶松弛前后的厚度变形结果
可将平均残余变形与试样原始厚度平均值的比值视为等效损伤塑性应变，丁腈橡胶和氢化丁腈橡胶不同压缩位移条件下的等效塑性应变见表3。
表3丁腈橡胶和氢化丁腈橡胶等效塑性应变
(3)辨识时变松弛-损伤本构模型的损伤塑性材料参数(k，r和)
通过在isight优化软件中集成abauqs搭建橡胶样件本构参数识别的优化模型。以时变松弛-损伤本构模型参数为优化设计变量，以试验数据和仿真数据之间的数值差异最小作为优化目标，优化模型的数学形式如下所示：
应用多岛遗传算法(miga)，以橡胶材料在不同压缩位移下的应力松弛试验数据为依据，识别本构模型的关键参数。多岛遗传算法是一种改进的并行分布遗传算法，与传统的遗传算法相比，它具有更好的全局求解能力、计算效率和收敛能力。结合优化模型，miga用于识别时变应力松弛-损伤本构模型参数的流程图如图7所示。具体识别流程如下：
首先，在寻优范围内随机生成整体种群，初始化之后将该种群分成若干个子种群，这些子种群在地理上分布在不同的岛内，子种群中个体携带的遗传信息中包含松弛参数k、r和不同阶段的等效应力每一次迭代过程均可融合实时的本构模型并进行有限元分析，提取并计算分析结果中的目标变量值；其次，计算每个子种群内个体的适应度值，按照迁移规则选取可迁移个体，实现与各个种群的联姻，可保证遗传信息的多样性；再次，各子种群内个体则按照一定的概率进行突变和交叉，如此循环往复，不停繁衍，不断更新优秀子种群，进而将优良的遗传信息写入新的橡胶材料本构模型中，直至仿真结果满足精度条件。
丁腈橡胶和氢化丁腈橡胶时变松弛-损伤本构模型的参数识别结果见表4。
表4丁腈橡胶和氢化丁腈橡胶损伤塑性材料参数(k，r和)
第五步：橡胶材料的时变松弛行为及损伤变形数值模拟预测
对比结果显示：识别出的本构模型参数能够准确描述橡胶的应力松弛退化行为及损伤变形效应，仿真结果与试验结果误差均小于5％，这说明数值模拟的预测结果准确性较高。
对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的，本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下，在其它实施例中实现。因此，本发明将不会被限制于本文所示的这些实施例，而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。