接触仿真知多少？详解ANSYS接触分析原理(一)
接触是一种很普遍的物理现象，它是Workbench用户最频繁使用的非线性特征之一。笔者在接触分析怎么搞？来看看这个ANSYS Workbench接触分析经典案例一文中，介绍了一个简单的案例。那么小伙伴们是否对ANSYS软件的接触解决方案的核心原理有所了解呢？
笔者就以本文整理一下接触分析的原理。本文篇幅较长，希望小伙伴们能有所收获，如有不当之处还请评论区指正~，也欢迎关注和支持本公 众号~
认识接触，接触有什么特征?
两个独立表面相互触碰并相切，称为接触。接触属于一种非线性问题，我们知道非线性问题主要有三种：几何非线性、状态非线性与材料非线性。而接触就是属于状态非线性问题。
当两接触体互相触碰或分离时，会发生刚度的突然变化，也就是状态改变的非线性行为。从现实物理现象上来讲，接触的这种状态变化具有以下特征:
① 不会渗透或穿透；
② 可传递法向压力和切向摩擦力；
③ 通常不传递法向拉力，在拉力下可自由分离。
接触的这些特征也是其仿真分析的变形协调要求。由于接触体间不相互渗透，所以程序必须建立两表面间的相互关系以阻止分析中的穿透。程序阻止穿透的行为我们称为强制接触协调。
此外，接触的状态(触碰、分离等)在外部边界条件下可能会不断变化，由于接触表面的法向刚度和切向刚度取决于接触或接触分离状态，这将导致接触刚度的不断变化，分析需跟踪接触状态，并对接触刚度进行更新，这会导致计算量加大，收敛困难。
因此，为了进行实为有效的计算，理解接触问题的特性和建立合理的模型就变得很重要。在利用软件进行接触分析时，其分析难点主要包括两个方面:
① 其一是分析过程中，我们并不知道接触区域是怎么变化的。因为接触区域受材料，载荷，边界条件等影响，其变化可能是瞬间变化的。
② 其二是大多的接触问题需要计算摩擦、接触压力，虽然ANSYS软件有几种摩擦等接触行为(无摩擦、粗糙、摩擦)和算法可供挑选，但它们也都是非线性的。计算分析时，除了需满足力收敛准则，位移收敛准则外，还需要满足接触收敛准则，因此摩擦、接触压力的计算会使问题的收敛性变得更复杂和困难，计算量也加大。
接触问题有哪些类型?
在接触问题分析中，分析模型中接触体可以是刚性的也可以是柔性的，刚性或者说刚度一般可以通过材料的弹性模量得到体现。通常，接触问题可按接触体的刚柔分类，并分为两类接触模型:
① 刚体一柔体接触
当一个表面除了允许刚体运动外，没有任何的应力、应变、和结构变形，这时我们可以认为该表面是一个理想的刚体表面。
但实际上这种完全刚性的表面是不存在的。相对的，当一个接触表面的刚度明显大于其它表面，且我们对该面的应力不感兴趣时，我们可以近似的将其作为刚体表面。
以上这种某一接触体可作为刚体的情况，就是刚体--柔体接触类型。
② 柔体一柔体接触
实际工程中，最常见的往往是两接触体都由柔性有限元单元构成的情况。如钢材与钢材接触，因为不同钢材的刚度差异小。这种情况，两个柔体的接触需要根据实际情况进行合理设置。
接触问题除按照以上刚柔分类外，还可按照接触方式划分，如:点——点接触，点——面(线)接触，面(线)——面(线)接触。这里线可以理解为二维的面。
接触体的刚柔、接触方式会影响接触算法对接触状态的检测，及目标面和接触面的选择。这些内容会在本文后续内容加以说明。
通常当两个任意形状的面接触时，应使用面一面接触单元，因为面-面接触具有以下特点:① 事先并不需要知道确切的接触位置；② 两接触表面可以使用不同的网格；③ 允许较大的相对滑动；④ 支持大应变和大转动支持材料非线性。这些特点可更好地模拟接触状态，有助于问题的收敛。
ANSYS如何处理接触分析，接触分析有什么流程?
在ANSYS 中，接触问题采用了接触对的概念进行处理，接触对由目标面和接触面组成。面上覆盖接触单元，就像皮肤一样铺设在有限元模型上。
接触面和目标面使用不同的接触单元类型。接触对通过实常数来识别，如下图所示。
图1 接触对的接触单元与目标单元
因此，不同的接触对必须定义不同的实常数组，即使实常数没有变化，也需要定义不同组号，因为如果一个问题中有多个接触区域，如果不定义不同的实常数就无法识别出接触区与和目标区的配对关系，接触定义就会发生混乱。
因此有几个接触对，就需要定义几个接触单元和目标单元以及其实常数组，即使所采用的接触单元一致。
那么ansys中接触单元的作用是什么呢？
由于接触表面的接触单元用于模拟真实的接触状态，故而必须满足接触的变形协调要求。因此，接触单元通常有以下三方面的作用:
① 防止接触面互相穿透(或使用接触协调)；
② 转换接触面之间的力传递(包括摩擦力和法向压力等)；
③ 对接触面的相对位置进行跟踪(接触位置、接触状态等)
为了实现上述接触单元的作用，在利用ANSYS进行接触分析时，首先需指定接触面和目标面定义接触，并建立和处理好初始的接触状态，然后根据实际情况设置好接触行为(或接触类型)，采用合适的接触数值算法计算。在ANSYS中，接触问题的分析过程通常如下:
图2 接触分析流程
04
如何定义接触面，如何进行接触检测？
前面提到接触分析前，需指定接触对的接触面和目标面，对于刚体-柔体接触，目标面总是刚体表面。对于柔体-柔体接触，接触面和目标面都看成是可变形的柔体。
由于接触算法中接触单元被限制不得穿透目标面(即接触面的接触单元积分点不能侵入目标面)，但目标面可以穿透接触面。如下图所示。
图3 接触单元积分点不能穿透目标面
这种限制使得柔-柔接触时，接触面和目标面的选择就需要遵循一定的准则:
① 如果凸面与平面或凹面接触，那么平面或凹面应该是目标面，凸面为接触面；
② 如果一个表面网格粗糙，而另一个表面网格较细，那么网格较粗的表面应该是目标面，细网格面为接触面；
③ 如果一个表面比另一个表面的刚度大(硬)，那么刚度大(硬)的表面应该是目标面，刚度小(软)的面为接触面；
④ 如果一个表面比另一个表面大，那么更大的表面应该是目标面，小面积的面为接触面；
⑤ 接触两侧分别为高阶单元对低阶单元，高阶单元面定义为接触面。
ANSYS在进行接触分析时，需判断接触面和目标面是否触碰，即接触处于何种状态，这就需要通过接触的跟踪和检测来实现。
接触的跟踪、检测需要通过接触检测点探明接触位置和接触状态。而接触检测点位于何处，可由接触单元的KEYOPT(4)选项定义，如下所示:
接触分析中，接触检测点是位于接触单元的积分点。积分点有节点积分和高斯积分两种(如下图4示)。接触单元的积分点不能侵入目标面,但是原则上目标面可以侵入接触面。
图5 节点积分点滑脱
但是有时必须采用基于节点的探测，例如尖角与线面的接触。如采用高斯积分点探测，节点与高斯积分点间会产生穿透，产生不准确的结果。
图6 节点积分
在ansys中，面-面接触单元直接使用接触面上的高斯积分点进行计算，而点-面接触单元直接使用接面上的节点进行计算(无高斯积分点)。
基于高斯积分点的探测默认指向接触面的法向；基于接触面节点的探测用于接触面比目标面光滑的情况；基于目标面节点的探测用于目标面比接触面光滑的情况。后两种基于节点的探测由于在探测前要计算接触面的法线方向，所以计算时间较基于高斯积分点的探测要长。
在ANSYS中，接触检测有“基于节点投影的接触（Nodal-Projection Normal from Contact）选项，如图7所示。它将接触面和目标面节点在法向投影的重叠区域强制定义接触约束，可产生如下效果:
① 对高阶单元结合 Normal Lagrange 法可以提供更精确的接触压力，且在接触边缘的接触压力和应变分布更加平滑；
② 对 Frictional 接触求解时可以很好地满足力矩平衡；但是不能与 MPC 接触匹配。
图7 接触面和目标面节点投影重叠区域
ANSYS如何确定接触状态？
ANSYS中，接触状态分为: 远离（Far）、接近（Near）、黏接（Sticking）和滑动（Sliding），如图8所示。
图8 接触状态
那么ANSYS 中，是如何确定接触状态的呢？在这里， 我们需要了解一下ANSYS 中的pinball区域(球形区域)。
因为pinball区域影响着接触分析的接触状态和一些其他接触参数的确定。根据ANSYS官方讲解，pinball区域的解释如下:
图9 ANSYS中的pinball区域的含义
pinball区域(弹球区域)区分的接触状态，为接触计算的效率提供了支撑，一般缺省设置对大多数问题有效，一般在软件中保持默认即可，但也可以根据需求修改。需要注意的是，程序只计算pinball区域内的接触穿透量，球形区域越大，程序所需的接触搜索时间越长，如下图10所示。
图16 接触穿透与接触间隙
图17 接触穿透
值得一提的是，真实状态下，物体在接触过程中是不允许穿透的，这是接触的基本特征，即接触变形协调所决定的。因此软件需要对穿透进行强制协调。
但在有限元分析过程中，如果不允许穿透，物体之 间发生接触或者取消接触时，会出现阶跃函数，导致收敛困难；为平滑处理接触状态的切换，如果允许一些极轻微的穿透， 接触将不再是一个突变函数，则较容易收敛，如图18所示。这种接触穿透的协调处理方式即软件中的罚刚度算法。
图18 接触状态过渡的处理
06
接触行为有哪些？对应哪些接触算法？
在ansys中，为提高接触的分析效率，软件根据接触的一些实际行为情况，提供了几种可供选择的接触表面行为。这些选项可以模拟许多不同的特殊物理效应。
这些接触行为包括: Bonded（绑定），No Separation（不分离），Frictionless（无摩擦），Rough（粗糙）和Frictional（摩擦）等。
具体分析时，这几种接触行为的表现如下:
绑定(Bonded)：法向不分离，切向无滑移。可以理解为目标面和接触面完全粘合在一起(默认设置），两个实体被认为焊接为了一个整体，是线性接触。如焊接部件一般为绑定接触。
不分离(No Separation)：允许接触面与目标面相对滑动，但是不允许接触面与目标面存在法向移动。即法向不分离，切向允许无摩擦的小滑移，也是线性接触。
无摩擦(Frictionless)：法向可分离，切向有无摩擦的滑移；是非线性接触。目标面和接触面可自由的分离和滑动。
粗糙(Rough):不 穿透、法向可分离、切向不滑移，是非线性接触；目标面和接触面间无滑移（和无穷大的摩擦系数类似）。
摩擦(Frictional)：允许有法向分离与切向有摩擦的滑移，这种接触类型应用的较多，但需要设置摩擦系数。
在进行仿真分析时，应根据不同的应用场合，判断是否存在法向分离和切向滑移，再根据切法向的情况选用合理的接触行为即可。
workbench中，不同接触行为有不同的数值计算算法作为分析的支撑，如表1所示。接触一般不允许穿透，为此定义了几种接触算法，以保证各种接触行为或状态的匹配，这等同于经典界面中的 keyopt（2）设置，即:
keyopt(2)=0，对应增广拉格朗日法
keyopt(2)=1，对应罚函数法
keyopt(2)=2，对应MPC法(多点约束)
keyopt(2)=3，对应法向拉格朗日，切向罚函数法
keyopt(2)=4，对应法向切向拉格朗日法
表1 不同接触行为的计算方法
不同的接触算法有不同的计算原理和特点，了解不同的算法原理，对正确选用合适的算法大有裨益。从图19可以看出，接触行为的算法分为三类五种，即罚刚度算法、多点约束算法、拉格朗日乘子法。其中，多点约束法(MPC)仅适用于绑定和不分离两种线性接触行为。
图19 接触算法分类及适用的接触行为
值得一提的是，ANSYS workbench中，拉格朗日乘子算法中的纯拉格朗日法(法向切向均为拉格朗日法)需通过插入命令流来实现。原因在于纯拉格朗日算法的收敛性不是很好，但一旦收敛，结果较为精确。
图20 workbench中的接触算法选择
在workbench中，通过插入command命令，利用单元KEYOPT关键字设置接触单元接触算法的方法如下图，VALUE值根据所需算法选择。这与在经典界面中的 keyopt（2）设置相通。
图21 workbench设置接触算法
07
接触算法是什么原理？
罚函数法是ANSYS 中的默认算法，适用于各类型的非线性接触行为(Frictional，Frictionless，Rough)，是相对于其他几种非线性算法中较为经济的一种算法。
罚函数法是将零件之间的接触假设成两个节点之间通过弹簧连接，通过以下计算公式来求解两个接触面之间的接触压力，罚函数(Pure Penalty)方程：
其中，Fn为法向接触力，Kn为法向接触刚度，Xp为穿透量。接触刚度Kn越高，穿透量Xp越小。如下图所示。
图22 罚函数法计算原理
理想情况下， Kn 为无限大，则 xp 为 0，但是如果Kn很大，会产生很大的接触反力，甚至让接触模型分离(颤振)，致使收敛困难或无法收敛；实际Kn取一个较大的数值， xp 较小以至忽略不计(经验上，一般取1E-8量级),也认为该方法可靠。这种方式实际上是通过改变罚刚度的值进行参数敏感性研究，从而对结果的有效性进行判断。
在实际情况下，两个零件表面是不会有穿透的，这是一种为增强收敛性而进行的数值近似方法，因此，穿透量越小，计算结果精度越高，但同时收敛性较差。因此，在使用罚函数算法的时候，需要仔细检查接触面的穿透量。
② 增强拉格朗日法（Augmented Lagrange）
增强拉格朗日算法是在罚函数的方法上衍生出来的一种方法，与罚函数法类似，但是在计算接触压力时，引入了附加项λ。增强拉格朗日方程：
增强拉格朗日方程因为有额外因子λ，使得增广拉格朗日法比罚函数法对接触刚度更不那么敏感。一般而言，接触刚度越高，穿透量越小。引入了λ之后，该算法下接触压力对于接触刚度的敏感性降低，更利于在给定的接触刚度较大的时候收敛，可以一定程度上提高计算精度，但是如果网格变形得过于扭曲，则计算迭代步数较多会造成收敛时间加长。
二、法向拉格朗日算法（Normal Lagrange）
法向拉格朗日算法中，是将接触压力作为一个自由度来满足接触兼容性，即Fn = 接触压力（DOF），如图23。因此它不需要计算接触刚度和穿透量来计算接触压力，而是将他看做一个自由度。
于是，有如下两种情况：零件不接触和零件接触。在计算过程中，这两种极限的情况会导致计算震荡剧烈从而较难收敛，但是一旦可以算收敛，由于这里没有假设零件之间的穿透，得到的结果精度较高。另外，拉格朗日法需要使用直接求解器来求解，计算速度较慢。
三、多点约束算法（MPC）
多点约束方式在内部增加了约束方程（ Constraint Equations）以绑定连接（Tie）接触面之间的位置，如图24所示。
图24 多点约束原理
这种方法直接有效地连接接触区域，而且可以适用于大变形开关开启的计算。主要用于 Bonded 和 No Separation。
特别适用于处理 Solid 与 Shell、 Shell 与 Shell 接触时易出现的接触面法向错误的情况。
例如，当软件提 示:“The normal of contact element XXX is not consistent with the normal of contact element XXX. Please use the ENORM command to correct it.”。
另外， MPC 算法是整个接触算法中求解速度最慢的，当 Solid 与 Shell 模型需要定义接触时，可以采用下面方法提高计算速度：①在 Solid 的表面手动新建一层非常薄的 Shell 模型；②Solid 模型与新建的 Shell 模型采用 Bonded 连接，新建的 Shell 模型与原来的 Shell 模型采用 Mesh Connections 连接。
四、接触算法的对比
对于绑定和不分离线性接触行为，由于是线性计算，其收敛性都比较好，计算速度也较快。对于三种非线性接触行为的计算，一般情况下，各算法从计算精度和收敛性上的排序可参考如下：
收敛性： 罚函数>增广拉格朗日>法向拉格朗日
精 度：一般拉格朗日>增广拉格朗日=罚函数
时 间：一般拉格朗日>增广拉格朗日>罚函数
对于个别情况，可能需要根据实际情况进行测试对比。各接触算法的特点如图25所示：
图25 接触算法的计算特征
可以看到，在接触的探测方面，纯罚函数和增强拉格朗日法默认基于高斯积分点的探测（On Gauss Points），一般较节点的探测更准确；拉格朗日和 MPC 法默认基于节点的探测（On Nodes-Normal from Contact 和 On Nodes-Normal to Target），较高斯积分点的探测点要少。
五、关于罚刚度算法中的刚度问题
罚刚度算法中，对增强拉格朗日法和罚函数法，需要法向和切向接触刚度。接触面和目标面之间的穿透量取决于法向刚度。在粘结接触中的滑动量取决于切向刚度。当两个物体接触时，接触刚度将被激活。
较高的刚度值可减少穿透量或滑动量，但会导致总体刚度知阵的病态和收敛困难。较低的刚度值会导致一定量的穿透/滑动，进而产生不准确的求解。理想地，使用足够高的刚度以便穿透/滑动量是少量可接受的，但是足够低的刚度在收敛方面是很好的。
在ANSYS中，法向接触刚度系数FKN、切向刚度系数FKT、穿透公差系数FTOLN、允许弹性滑移SLTO等实常数有缺省值，如图26所示。在大多数情况下，不需要定义这些接触刚度。
图26 接触刚度和接触穿透量的命令流实常数
在ANSYS workbench 中，可插入 Command，使用RMODIF设置接触实常数，命令如下：
RMODIF, NSET, STLOC, VALUE
其中:
NSET 表示接触对的实常数号，特定面可以用系统默认的CID 或TID；
STLOC 表示实常数序号，即上表中实常数表格的位置No.，对应具体实常数Name；其他接触实常数可查阅ANSYS help的接触实常数表(文末给出查询方法)；
VALUE 表示实常数定义的值。
如：RMODIF, CID, 3, 5 !表示将接触的法相接触刚度系数FKN改为5
在接触刚度控制方面，ANSYS Workbench在接触刚度控制中提供了两个图形界面选项，如下图27所示。
法向刚度如上图可在软件中直接设置，默认切向刚度是法向刚度的0.1，不能直接设置。法向刚度不能超过1E16，太大会因为计算机的舍入误差，导致求解精度下降。同时，计算时由于刚度太大会导致求解困难。一般实体之间的接触，接触刚度取默认值即可。但对于薄壁壳体，弯曲占主导，如冲压或滚压成型，刚度系数取值一般在0.01~0.1。
除法向刚度外，workbench中接触切向刚度不能直接设置。ANSYS自动定义了一个缺省的切向接触刚度，它与摩擦系数MU和法向刚度FKN成正比。默认切向刚度是法向刚度的10%，默认的切向刚度对应于默认值的切向刚度系数FKT=1.0。正的FKT值是因子，负的FKT值是切向刚度的绝对值。
ANSYS中，接触计算时接触刚度可随着计算更新，即ANSYS中的KEYOPT(10)=1或2，或KEYOPT(2)=3（法向拉格朗日乘子法和切向罚函数法）时，如下图28所示，ANSYS基于当前接触法向压力PRES和最大的允许弹性滑动SLTO，更新切向刚度FKT，即存在如下关系（根据节点计算出的刚度要乘上FKT）：
切向刚度FKT=×MU×PRES÷SLTO
=摩擦系数×接触法向压力÷最大的允许弹性滑动
当切向刚度FKT在每个迭代更新时，最大的允许弹性滑动SLTO实常数用于控制最大的滑动距离。ANSYS提供缺省的SLTO容差值，该值在大多数情况工作良好。但也可以覆盖SLTO的缺省值(接触对平均接触长度的1%)，采取自定义，但较大值将增强收敛但损害准确性。
图28 KEYOPT(10)和KEYOPT(2)选项
接触刚度对程序收敛和求解精度的影响最大。较大的刚度可提高求解精度，但使得收敛更加困难，因此必须谨慎定义接触刚度的大小。最适合的值和具体问题有关:
① 缺省值适用于大多数接触问题，但是，某些情形下程序提供的缺省值并不适合；
② 有时可能需要进行一些试验来取得一个既收敛又能保证求解精度的刚度值。
为确定一个合适刚度值，通常面临以下挑战：
① 使穿透最小以保证求解精度，因此，接触刚度应该非常大；
② 刚度太大将会有收敛困难，模型在接触表面可能来回振动，即出现接触颤振，如下图。
图29 接触在迭代计算时颤振示意
那么我们该如何在分析中，确定接触刚度和控制刚度计算呢？
通常在确定接触刚度时，确定合适的法向刚度，一般需改变法向刚度进行试算，查看穿透量及接触压力变化情况，直到接触压力稳定，穿透量复合容差要求，如下图所示。
在一个快速的大致的检查中，把模型的显示调为真实尺寸(true scale)后， 如果能观察到穿透现象，那么穿透可能过度了，这时应该增加刚度并重新开始分析； 如果需要过多的迭代或者根本不收敛，那么需要减小刚度值并重新开始分析； 重新分析时，使用接触刚度校正选项KEYOPT(10)来调整更新刚度值，这样可以在求解精度和收敛性之间取得一个较好的平衡。
好啦~，以上就是今天分享的全部内容了。关于软件接触分析的知识，始终需要回到软件本身，软件的帮助文档是最好的学习资料。本文对接触分析的理解也只能起到抛转引玉的作用，希望对小伙伴有所帮助，fighting~
