首页 > 模态试验 > 正文

有限元分析中的接触和摩擦模拟(三)

2018-10-11 21:04:50·  来源:模态空间  作者:北汽新能源王朋波  
 
摩擦数值算法简介在过去的三十年间,计算接触力学领域发展了多种用于求解摩擦接触问题的算法。在对摩擦接触问题进行数值模拟时,一个主要困难是摩擦力与切向滑动

 

摩擦数值算法简介

在过去的三十年间,计算接触力学领域发展了多种用于求解摩擦接触问题的算法。在对摩擦接触问题进行数值模拟时,一个主要困难是摩擦力与切向滑动量之间的本构关系是非光滑的,本构函数在某些点上不可微分,从而造成数值计算中迭代收敛困难,这个问题可以通过对摩擦本构关系的规则化来解决。

目前已有多种迭代方案用于带摩擦的接触分析,可大致分为以下几类:试探-校核算法、基于塑性理论中的弹塑性类比方法、基于优化理论的数学规划方法。后两类方法以严密的数学理论为基础,其可靠性高于试探-校核算法。

试探-校核算法通常是应用于小变形情况。对于摩擦接触问题,解的唯一性和存在性均不能保证,因此试探-校核算法在很多情况下并不可靠。但试探-校核算法仍然得到了成功的应用,在显式有限元分析中能够获得合理的结果。

近几十年,Coulomb定律和其他摩擦本构关系被纳入塑性理论的研究范畴,基于弹塑性理论的返回映射方案已成功应用于有限元摩擦接触分析,该方案使算法的收敛行为和可靠性产生了本质的提高。返回映射方案最初用于材料非线性问题,用以积分弹塑性本构关系。将摩擦定律类比为弹塑性本构关系,就可以将返回映射方案应用于带摩擦的接触分析。由于摩擦本构关系的非关联性,返回映射方案所得到的切线刚度矩阵通常不对称,增加了数值求解的难度。除返回映射方案之外,其他几种来源于塑性理论的方法,例如屈服极限拉氏乘子法,也已用于摩擦接触问题的数值分析中。

数学规划方案在摩擦的模拟中也有较为成功的应用。此外,其他的一些途径,例如内点法、屏障法、连续近似法等,也被尝试用于摩擦接触分析中。

摩擦接触分析中,需要区分滑动或粘结两种状态。对于粘结状态,通常可以使用一个参数(罚值)εT将其引入到罚函数或者拉格朗日格式,εT相当于切向滑动的刚度。但是这种处理使得静摩擦力的大小依赖εT的数值。对于过小的εT,可能导致物理上有问题的计算结果,即接触界面出现类似粘滑效应的行为,与实际情况不符。虽然可以采用较大的εT来解决这个问题,但是εT取值过大将导致方程组出现病态或者迭代流程失去健壮性。

8

摩擦定律与摩擦定律的正则化

 8.1 库伦摩擦定律

两接触表面的切向力一旦超过了某一个极限值,两表面不再保持粘结状态,而是发生相对的滑动。需要使用适当摩擦模型来描述这种切向行为。在实际工程中,库伦定律因其简单和适用性被广泛采用。库仑定律的表达式如下:

其中μ为摩擦系数,在经典库仑定律中被视为常数,由两接触面的材料种类决定;pN为法向接触压力;ġT为相对滑动速度;tT为切向摩擦力。

库仑定律表明,当摩擦力小于μpN时,无相对滑动,处于粘结状态;当有相对滑动时,摩擦力必然达到极限值μpN。

摩擦系数通常受多种因素的影响,例如表面粗糙度、相对滑动速度ġT、法向接触压力pN、环境温度θ等。如果考虑这些影响,可以得到经典库仑定律的改进形式:

 8.2 库仑定律的规则化

经典库伦摩擦定律是不可微的,在相对滑动速度ġT=0处发生由粘接至滑动的突变,当ġT反向时,摩擦力tT也立即发生反转。这种突变将导致数值计算中迭代的收敛困难。因此提出了一些规则化(即光滑化)的模型来代替经典库伦摩擦模型,这些规则化模型可以表达为如下形式

式中,函数φ(ġT)用来描述由粘结状态至滑动状态的光滑过渡。不同的规则化模型中,φ(ġT)可采用多种形式,分别对应于不同的规则化模型。

  • 平方根规则化模型

  • 双曲正切规则化模型

  • 反正切规则化模型
  • 分段直线模型

在以上模型中,标量参数ε为规则化参数,当ε→0即为经典库伦摩擦定律。ε不宜取值过大,否则规则化模型将与实际的接触界面粘结-滑动规律有较大出入。

上述规则化模型的前三种能够提供光滑的摩擦力-滑动量曲线,但第四种模型类似于理想弹塑性应力-应变关系,便于采用弹塑性力学中的相关方案处理。

当采用规则化的模型后,分析中不存在粘结状态,而统一的应用上面的数学描述来进行摩擦分析。这在物理上可以更真实地描述实际摩擦现象,因为两个接触面上都存在一定的不平度,因此完全粘结的接触状态是不存在的。只要有切向摩擦力的存在,总要伴随发生一定的相对滑动,而规则化的摩擦模型正好可以描述此类物理现象。

 8.3 摩擦定律的弹塑性类比

弹塑性类比方案的关键思想在于将切向滑动gT分离为弹性滑动部分gTe和塑性滑动部分gTp,如图2。

图2:切向滑移量的分离

这种分离方可视为对摩擦定律的规则化。但是它也有实际的物理含义:gTe代表接触界面上的弹性微小位移,卸载后gTe和由gTe引起的变形将恢复为0;gTp代表不可恢复的位移,在卸载后仍然存在。

对于切向摩擦的弹性部分,可采用以下各向同性线弹性模型,

式中,cT为界面弹性常数。

也可采用各向异性线弹性模型

式中,cT为弹性系数张量。

切向摩擦的塑性部分也有类似的本构方程,可根据弹塑性理论的标准概念导出。与弹塑性理论进行类比,可以建立屈服函数fs(tT)和相应的屈服准则,如下:

屈服函数fs(tT)可以简单的表达为:

上式等价于经典库伦摩擦定律。

利用最大耗散原理,得

然后可以导出塑性滑动的流动法则

此外,还可以得到加载-卸载准则

根据上式可确定参数ƛ.

9

商用有限元软件中的接触和摩擦算法

 9.1 ABAQUS中的接触算法

ABAQUS缺省使用拉氏乘子法施加接触约束,用户也可选择使用罚函数或增广拉格朗日法。在ABAQUS接触模拟中,要在各个构件上建立表面,定义可能会发生接触的一对表面(接触对),然后按照本构关系求出表面间的法向和切向相互作用。

ABAQUS可以自动确定二维和三维实体单元的自由面,并利用它们建立一个表面。对于壳单元和刚性单元,必须指明是正法向的面还是负法向的面。

ABAQUS应用单纯的主控-从属接触算法,执行单边接触检查:从属表面的节点不能穿透主控表面的任何部分,这种算法对于主控表面没有限制,它可以在从属表面的节点之间穿透从属表面。因此必须正确选择从属与主控表面,要遵循的简单原则是:

(1)从属表面应该是网格划分的更精细的表面。

(2)如果网格密度相近,从属表面应当由更柔软的材料组成。

定义接触对时,须确定两表面之间的相对滑动是小量还是有限值。默认设置是有限滑动状态。当应用小滑动公式时,从模拟开始就建立从属表面和主控表面之间的关系,确定主控表面哪个部分与从属表面的结点发生作用。这种相互关系在分析过程中保持不变。对于有限滑动情况,两个变形表面间的有限滑动模拟仅能用于二维问题,变形面和刚性面的有限滑动计算可用于三维问题。

图3为ABAQUS接触分析的流程。ABAQUS在增量步开始前检查所有接触对的状态,判断从属结点是否进入接触,对于接触的结点,判断是滑动还是粘接,然后施加相应的约束。

图3:ABAQUS接触算法流程

在检查力或力矩的平衡前,首先检查从属结点接触条件的变化,如果检测到当前迭代步的任何接触状态改变,则标识为严重不连续迭代,不进行平衡检验,而是改变接触约束进入第二次迭代。直至接触状态不再变化,才进行平衡收敛检查,如果收敛性不满足,则进行下一次迭代。

 9.2 MSC.MARC中的接触算法

MARC软件提供了基于直接约束的接触算法,是解决所有接触问题的通用方法。特别是对大面积接触,以及事先无法预知接触发生区域的接触问题,程序能根据物体的运动约束和相互作用自动探测接触区域。当发生接触时,使用边界条件直接约束运动体,两者的运动约束转化成了节点自由度的约束和节点力的约束。

MARC中的接触检查默认采用双边检查,为提高检测速度,用户也可选择单边检查。MARC不仅可以处理变形面和刚性面之间的有限滑动,对于变形面与变形面间的有限滑动也能很好的处理。

MARC软件接触算法基本流程如图4所示。从图4可见,整个过程包括:

  • 定义接触体。
  • 探测接触。
  • 施加接触约束。
  • 模拟摩擦。
  • 修改接触约束。
  • 检查约束的变化。
  • 判断分离和穿透。

此外,对于接触区域事先已知的情况,用户还可使用间隙单元(gap单元)方案,或者基于罚函数方法,通过用户子程序USPRNG对接触界面施加非线性弹簧。

图4:MARC接触算法流程

作者简介

王朋波,清华大学力学博士,汽车结构CAE分析专家。重庆市科协成员、《计算机辅助工程》期刊审稿人、交通运输部项目评审专家。专业领域为整车疲劳耐久/NVH/碰撞安全性能开发与仿真计算,车体结构优化与轻量化,CAE分析流程自动化等。

 
分享到:
 
反对 0 举报 0 收藏 0 评论 0
沪ICP备11026620号