非线性有限元本构模型(1)：非线性弹性模型

weixin

　　为了进行分析，选择材料模型是很重要，往往又不是很明确，仅有的信息可能是一般性的知识和经验，即可能是材料行为的几条应力-应变曲线。

　　在有限元软件库中选择合适的本构模型，如果没有合适的本构模型，要开发用户材料子程序。重要的是理解本构模型的关键特征，创建模型的假设，材料、荷载和变形域、以及程序中的数值问题是否适合模型。

　　应力-应变曲线
　　材料应力-应变行为的许多基本特征可以从一维应力状态(单轴应力或者剪切)的一组应力-应变曲线中获得，多轴状态的本构方程常常基于在试验中观察到的一维行为而简单生成。

　　定义伸长

　　名义应力(工程应力)给出为：

　　工程应变定义为：

　　载荷-位移曲线

　　以每单位当前长度应变的增量随长度的变化得到另一种应变度量。

　　对数应变(也称为真实应变)

　　对材料时间求导，表达式为：

　　在一维情况下，上式为变形率。

　　当前面积的表达式给出为：

　　Cauchy(或者真实)应力表示为：

　　工程应力-应变曲线

　　真实应力-应变曲线

　　考虑一种不可压缩材料(J=1)，名义应力和工程应变的关系为：

　　真实应力(对于不可压缩材料)

　　说明了对于本构行为应用不同泛函表达式的区别，对于同样材料取决于采用何种应力和变形的度量。

　　应力-应变曲线的显著特征之一是非线性的度。材料线弹性行为的范围小于应变的百分之几，就可以采用小应变理论描述。

　　应力-应变反应与变形率无关的材料称为率无关;否则，称为率相关。名义应变率定义为：

　　因为：

　　即名义应变率等于伸长率，例如：

　　可以看出，对于率无关材料的应力-应变曲线是应变率独立的，而对于率相关材料的应力-应变曲线，当应变率提高时是上升的;而当温度升高时是下降的。

　　率无关和率相关材料的一维反应

　　对于弹性材料，应力-应变的卸载曲线简单地沿加载曲线返回，直到完全卸载，材料返回到了它的初始未伸长状态。对于弹-塑性材料，卸载曲线区别于加载曲线，其斜率是典型的应力-应变弹性(初始)段的斜率，卸载后产生永久应变。其它材料的行为介于这两种极端之间。由于在加载过程中微裂纹的形成材料已经损伤，脆性材料的卸载行为是当荷载移去后微裂纹闭合，弹性应变得到恢复。卸载曲线的初始斜率给出形成微裂纹损伤程度的信息。

　　弹性

　　弹-塑性

　　弹性含损伤

　　一维弹性
　　弹性材料的基本性能是应力仅依赖于应变的当前水平。这意味着加载和卸载的应力-应变曲线是一致的，当卸载结束时材料恢复到初始状态。称这种应变是可逆的。而且，弹性材料是率无关的(与应变率无关)。弹性材料的应力和应变是一一对应的。

　　小应变

　　对于任意应变，不管如何达到应变值，上式给出唯一应力值。

　　可逆和路径无关默认在变形中没有能量耗散，在弹性材料中，储存在物体中的能量全部消耗在变形中，卸载后材料恢复。

　　对于一维弹性材料，可逆、路径无关、无能量耗散是等价的特征。

　　对于二维和三维弹性，以及超弹性材料，也类似。

　　应变能一般是应变的凸函数，例如，

　　如有：

　　则公式的等号成立。

　　凸应变能函数的一个例子如图所示。在这种情况下，函数是单调递增的;如果w 是非凸函数，则s先增后减，材料应变软化，这是非稳定的材料反应，ds/dεx<0如右下图。

　　凸应变能函数 应力应变曲线

　　非凸应变能函数 相应的应力应变曲线

　　大应变
　　从弹性推广到大应变，只要选择应变度量和定义应力(功共轭)的弹性势能。势能的存在是默认了可逆、路径无关和无能量耗散。如：

　　在弹性应力-应变关系中，从应变的势函数可以获得应力为超弹性。如一维大应变问题，以Green应变的二次函数表示

　　对于小应变问题，即为胡克定律。

　　一种材料的Cauchy应力率与变形率相关，称为次弹性。这种关系一般是非线性的，给出为：

　　一个特殊的线性次弹性关系给出为

　　对上式的关系积分，得到

　　这是与路径无关的超弹性关系。对于多轴问题，一般次弹性关系不能转换到超弹性，它仅在一维情况下是严格路径无关的。然而，如果是弹性小应变，其行为足以接近路径无关的弹性行为。因为次弹性的简单性，次弹性公式的多轴一般形式常常应用在有限元软件中，以模拟大应变弹塑性的弹性反应。

　　非线性弹性
　　对于有限应变有许多不同的应力和变形度量，同样的本构关系可以写成几种不同的形式，总是可能从一种形式的本构关系转换到另一种形式。

　　大应变弹性本构模型首先表述成Kirchhoff材料的一种特殊形式，由线弹性直接生成到大变形。满足路径无关、可逆和无能量耗散。因此，路径无关的程度可以视为材料模型弹性的度量。

　　次弹性材料是路径无关程度最弱的材料，遵从Cauchy弹性，其应力是路径无关的，但是其能量不是路径无关的。

　　超弹性材料或者Green弹性，它是路径无关和完全可逆的，应力由应变势能导出。

　　小应变和大转动
　　许多工程应用包括小应变和大转动。在这些问题中，大变形的效果主要来自于大转动，如直升机旋翼、船上升降器或者钓鱼杆的弯曲。由线弹性定律的简单扩展即可以模拟材料的反应，但要以PK2应力代替其中的应力和以Green应变代替线性应变，这称为Saint-Venant- Kirchhoff材料，或者简称为Kirchhoff材料。最一般的Kirchhoff模型为

　　式中C为弹性模量(切线模量)的四阶张量，对Kirchhoff材料是常数，代表了应力和应变的多轴状态。它可以完全反映材料的各向异性。

　　一般的四阶张量有3∧4=81个独立常数，与全应力张量的9个分量和全应变张量的9个分量有关。 如次弹性本构方程：

　　式中C为弹性模量的四阶张量，有81个常数。利用对称性可以显著地减少常数。

　　利用势能表示的应力-应变关系和Green公式，

　　故有

　　这样C为对称矩阵(主对称性)，在81个常数中有45个是独立的。成为上三角或下三角矩阵。

　　应力张量和应变张量均为对称张量(次对称性)，即

　　应力和应变张量的对称性要求应力的6个独立分量仅与应变的6个独立分量有关，由弹性模量的局部对称结果，独立常数的数目减少到36个。

　　再利用模量的主对称性使独立弹性常数的数目减少，由36个常数减少为21个，为各向异性材料。

　　写成矩阵形式为(可以是上或下三角矩阵)

　　对于正交各向异性，具有正交的三个弹性对称面，当坐标变号，为使应变能密度不变，有

　　这样由21个常数减少为14个，为正交各向异性材料。若材料对称坐标平面，当沿轴平面反射时，弹性模量不变，因为正交各向异性体，有

　　对于一个由三个彼此正交的对称平面组成的正交材料(如木材或纤维增强的复合材料)，仅有9个独立弹性常数，Kirchhoff应力-应变关系为材料对称坐标平面，为正交各向异性体

　　对于各向同性材料，仅有3个常数

　　材料对称的一个重要的例子是各向同性。一个各向同性材料没有方位或者方向的选择，因此，当以任何直角坐标系表示的应力-应变关系是等同的。对于小应变的许多材料(如金属和陶瓷)可以作为各向同性进行模拟。张量C是各向同性的。在任何坐标系统中，一个各向同性张量有相同的分量。

　　对于各向同性的Kirchhoff材料，其应力-应变关系可以写成为

　　式中Lamé常数，体积模量K，杨氏模量E和泊松比v的关系为

　　不可压缩性
　　在变形的过程中，不可压缩材料的体积不变，密度保持常数。不可压缩材料的运动称为等体积运动。

　　总体变形为：

　　等体积约束运动的率形式：

　　将应力和应变率度量写成偏量和静水(体积的)部分的和，对于不可压缩材料，静水部分也称为张量的球形部分，分解式为：

　　对于不可压缩材料，压力不能从本构方程确定，而是从动量方程确定。

　　Kirchhoff应力

　　由Jacobian行列式放大，称它为权重Cauchy应力。对于等体积运动，它等同于Cauchy应力。

　　次弹性
　　次弹性材料规律联系应力率和变形率。

　　上式是率无关、线性增加和可逆的。对于有限变形状态的微小增量，应力和应变的增量是线性关系，当卸载后可以恢复。然而，对于大变形能量不一定必须守恒，并且在闭合变形轨迹上作的功不一定必须为零。次弹性规律主要用来代表在弹-塑性规律中的弹性反应，小变形弹性，且耗能效果也小。

　　切线模量之间的关系
　　某些次弹性本构关系共同应用的形式为

　　对于各向同性材料Jaumann率的切线模量为

　　对于同一种材料，切线模量不同，材料反应的率形式不同，如

　　对于简单剪切问题采用相同的材料常数使用不用应力客观率的应力比较

　　超弹性材料
　　平衡方程是以物体中应力的形式建立的，应力来源于变形，如应变。如果本构行为仅是变形的当前状态的函数，为与时间无关的弹性本构。而对于接近不可压缩的材料，仅依赖变形(应变)不一定能够得到应力。

　　储存在材料中的能量(功)仅取决于变形的初始和最终状态，并且是独立于变形(或荷载)路径，称这种弹性材料为超弹性(hyper-elastic)材料，或者为Green弹性，例如常用的工业橡胶。动物的肌肉也具有超弹性的力学性质。这里主要讨论橡胶材料的超弹性力学行为。

　　对于功独立于荷载路径的弹性材料称之为超弹性(Green弹性)材料。超弹性材料的特征是存在一个潜在(或应变)能量函数，它是应力的势能：

　　通过适当转换获得了对于不同应力度量的表达式

　　由于变形梯度张量F是不对称的，因此名义应力张量P的9个分量是不对称的。

　　在橡胶大变形中应用多项式模型和Ogden指数模型。

　　橡胶是提取橡胶树、橡胶草等植物的胶乳，加工后制成的具有弹性、绝缘性、不透水和空气的材料。在半个世纪前，“橡胶”一词是专指生橡胶，它是从热带植物巴西三叶胶的胶乳提炼出来的。

　　目前，世界半数以上的橡胶是合成橡胶。合成橡胶的种类很多，例如，制造轮胎使用的丁苯橡胶(苯乙烯和丁二烯的共聚物)或乙丙烯橡胶(ERP);用于汽车配件的有氯丁橡胶及另一种具有天然橡胶各种性能的异戊橡胶。

　　在众多的合成橡胶中，硅橡胶是其中的佼佼者。它具有无味无毒，不怕高温和严寒的特点，在摄氏300度和零下90度时能够“泰然自若”、“面不改色”，仍不失原有的强度和弹性。例如生物材料。

　　橡胶具有许多特殊的性能，例如电绝缘性、耐氧老化性、耐光老化性、防霉性、化学稳定性等。

　　1839年，Charle Goodyear发明了橡胶的硫化方法，其姓氏现在已经成为国际上著名橡胶轮胎的商标。

　　从19世纪中叶起橡胶就成为一种重要的工程材料。然而，橡胶材料的行为复杂，不同于金属材料仅需要几个参数就可以描述材料特性。橡胶材料受力以后，变形是伴随着大位移和大应变，其本构关系是非线性的，并且在变形过程中体积几乎保持不变。

　　橡胶是一种弹性聚合物，其特点是有很强的非线性粘弹性行为。它的力学行为对温度、环境、应变历史、加载速率都非常敏感，这样使得描述橡胶的行为变得非常复杂。橡胶的制造工艺和成分也对橡胶的力学性能有着显著的影响。

　　由于计算机以及有限元数值分析的飞速发展，我们可以借助计算机来对超弹性材料的工程应用进行深入研究以及优化设计。可以用有限元等数值方法来计算分析橡胶元件的力学性能，包括选取和拟合橡胶的本构模型，以及用有限元建模和处理计算结果等。

　　常用的橡胶性态可分为固体橡胶和泡沫橡胶。固体橡胶是几乎不可压缩的，其泊松比接近于0.5。可逆，大应变。初始各向同性，应变增加后分子定向排列。

　　固体橡胶材料的拉伸试验曲线与材料演化模型

　　一般将多孔橡胶或弹性泡沫材料统称为泡沫材料。弹性泡沫材料的普通例子有多孔聚合物，如海绵、包装材料等。

　　泡沫橡胶是由橡胶制成的弹性泡沫材料，能够满足非常大的弹性应变要求，拉伸时的应变可以达到500%或更大，压缩时的应变可以达到90%或更小。与固体橡胶的几乎不可压缩性相比，泡沫材料的多孔性则允许非常大的体积缩小变形，因此具有良好的能量吸收性。

　　泡沫橡胶材料的多面体微元模型

　　a) 开放腔室，b) 封闭腔室

　　小应变 <5%，线弹性，泊松比为0.3 。大应变，压缩时，泊松比为0.0;拉伸时，泊松比大于0.0。

　　泡沫橡胶材料的应力-应变曲线

　　a)压缩 b)拉伸

　　橡胶本构模型

　　典型固体橡胶材料单轴拉伸应力-应变曲线

　　小变形
　　以多项式形式本构模型为例，其应变能密度表达式为：

　　忽略二阶及二阶以上小量，变为：

　　弹性常数为：

　　当：

　　定义伸长

　　工程应变定义为：

　　二阶张量基本不变量

　　小变形，有：

　　常用的橡胶力学性能描述方法主要分为两类，一类是基于热力学统计的方法，另一类是基于橡胶为连续介质的唯象学描述方法。

　　热力学统计方法的基础为观察到橡胶中的弹性恢复力主要来自熵的减少。橡胶在承受荷载时分子结构无序，熵的减少是由于橡胶伸长使得橡胶结构由高度无序变得有序。由对橡胶中分子链的长度、方向以及结构的统计得到本构关系。

　　唯象学描述方法假设在未变形状态下橡胶为各向同性材料，即长分子链方向在橡胶中是随机分布的。这种各向同性的假设是用单位体积(弹性)应变能函数(U)来描述橡胶特性的基础，其本构模型为多项式形式模型和Ogden形式模型。

　　典型的本构模型为多项式形式，其应变能密度表达式为：

　　特殊形式可以由设定某些参数为0来得到。如果所有

　　则得到减缩多项式模型

　　对于完全多项式，如果N =1 , 则只有线性部分的应变能量，即Mooney-Rivlin形式

　　对于减缩多项式，如果 N =1 ，则得到Neo-Hookean形式

　　Mooney-Rivlin形式和Neo-Hooken形式本构模型

　　(后者是将Hooke定律扩展至大变形)

　　Yeoh形式本构模型是 N = 3 时减缩多项式的特殊形式

　　典型的S形橡胶应力-应变曲线 ，C10正值，在小变形时为切线模量;C20为负值，中等变形时软化;C30正值，大变形时硬化。

　　其他形式的本构模型有：
　　Ogden形式本构模型

　　Arruda-Boyce形式本构模型

　　Van der Waals模型

　　试验拟合本构模型系数：

　　橡胶类材料的本构关系除具有超弹性、大变形的特征外，其本构关系与生产加工过程有直接关系，如橡胶配方和硫化工艺。确定每一批新加工出来的橡胶的本构关系，都要依赖于精确和充分的橡胶试验。

　　通常在试验中应该测得在几种不同荷载模式下的应力-应变曲线，这样可以选择出最合适的本构模型以及描述这种模型的参数。

　　同一种橡胶材料的三种拉伸变形状态的应力-应变曲线图，对比试验曲线，由最小二乘法拟合多项式本构模型中的系数。

　　给出实验数据，应力表达式的系数通过最小二乘法拟合确定，这样可以使得误差最小。即对于n 组应力-应变的试验数据，取相对误差E 的最小值，拟合应力表达式中的系数，得到理论本构模型。

　　确定材料常数的经验公式

　　对于已经成型的橡胶元件，通常不容易通过上述试验来确定其材料常数。经验公式是通过橡胶的IRHD硬度指标来确定材料的弹性模量和切变模量，再由材料常数和弹性模量的关系来确定材料常数。基本公式为(小应变条件)

　　将得到的材料常数代入Mooney-Rivlin模型进行计算。

　　结论与讨论
　　由于大型有限元软件的迅速发展，使得复杂的超弹性模型计算过程由计算机程序完成，在ABAQUS等商用软件中给出了具体的计算。用户要熟悉如何输入数据文件，根据试验数据拟合和选用合适的本构模型，如何处理输出结果并检验其是否正确。对于初学者来说，商用软件是一个“黑匣子”，因此，掌握超弹性材料模型理论和计算方法是取得仿真成功的关键。

　　需要注意的是，对于不可压缩材料的平面问题，无论是解析解还是数值解，均不能采用平面应变解答。因为对于不可压缩材料，如果采用平面应变模型，其体积不变，内力为不确定量，在有限元中的节点位移不能反映单元内力的变化。对于不可压缩材料或者接近于不可压缩材料的平面问题，务必应用平面应力(或者广义平面应变)解答。

　　本文根据百度文库PPT讲义《非线性专题-本构模型》整理而成。讲义原作者庄茁，由qiannianyiri分享。