分享|接触分析

ANSYS有三种类型的接触单元：

点对点：最终位置事先知道；只能用于低次单元

点对面：接触区域未知，并且允许大滑动；

面对面：接触区域未知，并且允许大滑动（相对点对面接触有几个优点）。

接触分析属于高度非线性分析，需要较多的计算资源，这对网格划分以及接触面的选择提出了较高要求。

ansys可完成的接触分析主要有三类：点点，点面，面面接触分析；

接触分析主要分为两类：刚体—柔体接触以及柔体—柔体接触。其中，金属成型分析是典型的刚体柔体接触，一般的接触的问题均为柔体——柔体接触。
★分析的难点在于：1.接触面的识别和选择；2.摩擦模型的选择。
ansys接触分析是通过建立一层接触单元覆盖在接触面之上

点点接触一般较少使用，它适用于：预先知道接触位置，且相对滑动忽略，转动量很小，即使是几何非线性分析。一些过盈装配问题可以采用点点接触代替面面接触；
点面接触不需要知道确切接触位置，也不必保持网格一致，允许较大的变形和相对滑动。这种接触推荐采用contact48而不是26来计算；
面面接触是最为常见也是适用范围较广的接触类型：
★几个原则（asymmetric contact）：接触单元不能渗入目标面，但是目标（面上的）单元可以渗入接触面。目标面总是刚性的，接触面总是柔性的。
平面或者凹面为目标面；
网格细致的为接触面，网格粗糙的为目标面（目标面可以被渗入）；

The softer surface should be the contact surface and the stiffer surface should be the target surface.
高阶单元为接触面，低阶为目标面；However, for 3-D node-to-surface contact, 低阶单元为接触面，高阶为目标面；
面积大的是目标面。

In the case of 3-D internal beam-to-beam contact modeled by CONTA176 (a beam or pipe sliding inside another hollow beam or pipe),内部的为接触面，外部为目标面；However, when the inner beam is much stiffer than the outer beam, the inner beam can be the target surface.
若不能很好的区分接触面和目标面（When there are several contact pairs involved in the model, and the graphical picking of contact and target surfaces is difficult, you can just define the symmetric contact pairs and, by setting KEYOPT(8) = 2）可采用对称接触分析（Symmetric Contact），即通过设置KEYPOINT(8)=2 实现。但是这样效果并不会很理想，计算消耗更大。
★接触分析属非线性程度较高，宜采用位移加载方式，而不是力加载方式以确保收敛。事实上，施加力很难收敛；
接触面和目标可采用ansys自带接触管理器(Contact Manager)生成，非常方便；
→接触管理器只能识别的六面体网格（四边形网格），有一定的局限性。对于自由划分的实体不能使用接触管理器；
★因为大的时间增量会使因为大的时间增量会使因为大的时间增量会使因为大的时间增量会使迭代趋向于变得不稳定，使用线性搜索选项来使计算稳化。

打开线性搜索，打开自动时间步都有助于收敛，但不能与弧长法同时使用，否则在激活弧长法时发出警告并失效。LNSRCH,1;AUTOTS,ON；

非线性分析预测器，命令：PRED，Sskey，——，Lskey

该命令作用是加速收敛，在大转动或黏弹性分析中可能导致发散，需关闭。（当存在转角自由度或使用SOLID65单元时缺省为关闭）
To avoid a slow convergence rate and use an updated stiffness matrix,采用完全的牛顿—拉夫逊迭代方法。NROPT,FULL,,OFF
接触刚度宜在0.01~1之间取值，先取较小值。很多情况下，不收敛是由于选择了过大的接触刚度。

设置合理的平衡迭代次数，一个合理的平衡迭代次数通常在25和50之间。

命令：NEQIT GUI：Main Menu>Solution>-load step opts-Nonlinear>Equilibriwm iter                                             

 其中常用的设置有：

KEYPOINT(2)——接触算法；

接触刚度，FKN；
摩擦类型,；
KEYPOINT(4)——接触方位点的位置；
初始接触条件；
KEYPOINT(7)——自动时间预测控制；
KEYPOINT(9)——初始渗透的影响；
KEYPOINT(12)——接触表面情况。
KEYPOINT(2)——接触算法：a 采用默认的拉格朗日法是比较合理的，使用拉格朗日算法的同时应使用实常数FTOLN，默认为0.1，这是一个收敛准则不宜过大。
接触刚度，FKN：接触刚度的意义是，两个表面渗透量的大小。接触刚度越大则渗透越小。但是这个值不宜过大，会导致不收敛；也不能过小，会导致渗透过大结果失真，这个值在0.01到10之间，默认为1.0。一般来说，先选取一个较小的接触刚度，保证前几个子步能够收敛。结合KNF和FTOLN调整。
摩擦类型,：d 程序默认接触面之间无摩擦，默认取最大摩擦力限制为：TAUMAX=1.2e20 。实际上此值可取：根号3分之屈服应力。这个地方存疑。
KEYPOINT(4)——接触方位点的位置：s 有两类方式：高斯积分点和牛顿——cots积分点。其中高斯积分点比较准确，且第二种只对于焦点接触有效果。
初始接触条件：s 这一点很可能跟为什么力加载不容易收敛有关系！ 是在静力分析中，接触分析很容易出现刚体位移约束不足的情况，这也是为何加位移荷载易收敛的一个原因。因为在建模过程中很难保证两个接触面完全的贴合，接触不能有效的约束，导致发生刚体位移。解决办法尝试：1.设置ICONT来定义一个接触环，这个值不宜过大，接触环内的节点被强制接触在一起；2.设置是常数PMAX,PMIN 来指定初始容许渗透范围，这样当分析是渗透若果超越PMAX，程序会调整接触面减小渗透，是否可以考虑设定一个负的PMIN来强制接触节点 有待考证h或许可以解决接触问题加力的问题；3.设置KEYPOINT(9)：KEYPOINT(9)=0，包含初始渗透；KEYPOINT(9)=1，则消除初始渗透；KEYPOINT(9)=2，过度渗透逐渐过度要0，一般用于过盈分析中；
KEYPOINT(12)——接触表面情况：z代表接触面之间的分离类型。
0 -- Standard，单边法相接触，可分离
1 -- Rough，法相接触，但是无滑动，不可分离
2 -- No separation (sliding permitted)，不分离，可滑动
3 -- Bonded
4 -- No separation (always)
5 -- Bonded (always)
6 -- Bonded (initial contact)

 反复荷载下约束混凝土本构关系模型与应用  赵文艳; 张素梅  大庆石油学院学报 

利用ANSYS Solid65单元分析混凝土结构  霍越群; 赵光宇  吉林省土木建筑学会2011年学术年会论文集 

1.

tb,concr,2  

tbtemp,0

tbdata,,0.5,1,2.868,-1这几个的值比较重要,要仔细反复调整.

 2.计算中最好关闭压碎,不关闭也可以，不过得到的曲线不光滑.

 接触刚度问题汇总：

所有的接触问题都需要定义接触刚度，两个表面之间渗 量的大小取决了接触刚度，过大的接触刚度可能会引起总刚矩阵的病态，而造成收敛困难，一般来谘，应该选取足够大的接触刚度以保证接触渗透小到可以接受，但同时又应该让接触刚度足够小以使不会引起总刚矩阵的病态问题而保证收敛性。

   程序会根据变形体单元的材料特性来估计一个缺省的接触刚度值，你能够用实常数FKN来为接触刚度指定一个比例因子或指定一个真正的值，比例因子一般在0.01和10之间，对于大变形问题，通常为0.01-0.1，默认为1.0。用户应当总是检验以使穿透到极小值，而又避免过多的叠代次数。当避免过多的迭代次数时，应该尽量使渗透到达极小值。 FTLON的缺省值为0.1，你可以改变这个值但要注意如果此值太小可能会造成太多的迭代次数或者不收敛

为了取得一个较好的接触刚度值，又可需要一些经验，你可以按下面的步骤过行。

1、开始时取一个较低的值，低估些值要比高估些值好因为由一个较低的接触刚度导致的渗透问题要比  过高的接触刚度导致的收敛性困难，要容易解决。

2、对前几个子步进行计算，直到最终载荷的一个比例（刚好完全建立接触）。

3、检查渗透量和每一子步中的平衡迭代次数，如果总体收敛困难是由过大的渗透引起的（而不是由不平衡力和位移增量引起的），那么可能低估了FKN的值或者是将FTOLN的值取得大，如果总体的收敛困难是由于不平衡力和位移增量，达到收敛值需要过多的迭代次数，而不是由于过大的渗透量，那么FKN的值可能被高估。

4、按需要调查FKN或FTOLN的值，重新分析。

5、调整接触刚度同时调整接触区单元密度(或大小)才能更容易解决问题。从我的经验来看，合适的单元密度更能保证解的准确！因为很多接触问题仅靠调整接触刚度虽然能得到的解，但并不准确。

在有限元分析中，接触单元通常用来描述两物体相互接触或滑动的界面。近年来，ANSYS开发了一系列的接触单元。刚开始有节点对节点单元CONTAC12和CONTAC52，接着有节点对地单元CONTAC26，然后有节点对面单元CONTAC48和CONTAC49。最近几年，我们引入一类面对面接触单元CONTA169和CONTA174，同时还有一种新的节点对节点单元CONTA178。 
虽然接触单元的参数具有多样性，但我们在使用他们时可谨记重要的一点，他们具有一个共同的特点，即除了CONTA178的KEYOPT（2）=0或1外，所有的接触单元都有接触刚度。在现实中实际上相邻结构之间只是一种空隙，但在有限元分析中，这种空隙是一带有刚度的接触单元，这是因为通过刚度矩阵来实现接触算法的。一些接触单元要求使用者输入刚度值，同时另外的接触单元若没有输入则使用缺省值。分析工程师所面对的问题就是针对给定的条件确定一个合理的刚度值。如果过高，问题将会不收敛，如果过低，可能得到错误的结果。那么我们所面对的问题是怎样才能找到一个正确的刚度值？ 
我认为唯一的方法就是我们必须试用不同的值直到找到正确的值。也就是刚开始我们应该使用一个较小的值，然后稳步的增加直到分析的结果不再有什么变化。那么对于我们这一特定分析的问题，这一点就是我们所想要的合适值。 
我们可举例说明，如图1所示，平行放置两个悬臂梁，并有少许的交迭，下面的左边固支，上面的右边固支，当在上面梁的自由端施加一个向下位移时，梁变形弯曲并接触下面的梁，然后一起向下运动。用SOLID45单元划分梁，用TARGE170和CONTA174面面接触单元来描述相互作用。在此基础上，把CONTA174单元的刚度从非常低变到非常高，从而来观察它对结果的影响和收敛的迭代次数。图2说明了下梁自由端的偏移随接触单元刚度的变化情况，当刚度增加时，偏移量接近一个常数值（我们可以假定它是一个'正确'的结果。）图3说明求解所需的迭代次数，当接触单元刚度增加时，求解所需的迭代次数也是增加的，并服从指数关系。如果刚度过高，问题很有可能根本就不收敛。图4说明在上梁自由端接触单元的渗透量，当刚度增加时，渗透量降低。 
从这些图可知，当接触单元的刚度为10e6时，可获得合理精确的结果。任何大于该值的刚度对下梁的偏移量没有什么影响，而求解所需的迭代数却显著的增加。对于这个题目，10e6的刚度是很适合的。但是，如果改变边界条件、网格密度、两梁之间的相对位置、材料特性或梁的几何形状，能获得满意结果的接触刚度值将是不同的。比如，如果网格密度增加，则接触单元数将增加，每一个单元上的载荷将降低。如果接触单元数增加两倍，一个合适的接触单元刚度值应为原来的一半。 
由于每个题目都是不一样的，所以在求解之前并没有通用的方法来确定接触单元刚度的最佳值。我们不得不试算一个我们认为合适的值然后查看计算结果。一个有经验的分析工程师可能只查看一个计算结果来判定所取值的合适度，但对于大多数情况而言，最好用一个合理而不过度精确的刚度值进行第一次求解，然后用10倍于该值的刚度进行第二次求解，如果两者结果相差很小，而迭代数增加很多，那么我们则正好取得了曲线上的突变点，从而获得相当好的结果。 
接触单元刚度问题仅仅是一个例子，即对于分析工程师来说，总是置疑于分析结果的正确与否是非常重要的，并要意识到数值仿真的局限性和潜在的假设及他们怎样影响所分析问题的结果。 
图1

接触问题的计算方法

接触问题的关键在于接触体间的相互关系（废话 ），此关系又可分为在接触前后的法向关系与切向关系。 

法向关系： 

在法向，必须实现两点：1）接触力的传递。2）两接触面间没有穿透。 

ANSYS通过两种算法来实现此法向接触关系：罚函数法和拉格朗日乘子法。 

1．罚函数法 

是通过接触刚度在接触力与接触面间的穿透值（接触位移）间建立力与位移的线性关系： 

接触刚度*接触位移=法向接触力 

对面面接触单元17*，接触刚度由实常数FKN来定义。 穿透值在程序中通过分离的接触体上节点间的距离来计算。接触刚度越大，则穿透就越小，理论上在接触刚度为无穷大时，可以实现完全的接触状态，使穿透值等于零。但是显而易见，在程序计算中，接触刚度不可能为无穷大（否则病态），穿透也就不可能真实达到零，而只能是个接近于零的有限值。 

以上力与位移的接触关系可以很容易地合并入整个结构的平衡方程组K*X=F中去。并不改变总刚K的大小。这种罚函数法有以下几个问题必须解决： 

1）接触刚度FKN应该取多大？ 

2）接触刚度FKN取大些可以减少虚假穿透，但是会使刚度矩阵成为病态。 

3）既然与实际情况不符合的虚假穿透既然是不可避免的，那么可以允许有多大为合适？ 

因此，在ANSYS程序里，通常输入FKN实常数不是直接定义接触刚度的数值，而是接触体下单元刚度的一个因子，这使得用户可以方便地定义接触刚度了，一般FKN取0.1到1中间的值。当然，在需要时，也可以把接触刚度直接定义，FKN输入为负数，则程序将其值理解为直接输入的接触刚度值。 

对于接近病态的刚度阵，不要使用迭代求解器，例如PCG等。它们会需要更多的迭代次数，并有可能不收敛。可以使用直接法求解器，例如稀疏求解器等。这些求解器可以有效求解病态问题。 

穿透的大小影响结果的精度。用户可以用PLESOL,CONT,PENE来在后处理中查看穿透的数值大小。如果使用的是罚函数法求解接触问题，用户一般需要试用多个FKN值进行计算，可以先用一个较小的FKN值开始计算，例如0.1。因为较小的FKN有助于收敛，然后再逐步增加FKN值进行一系列计算，最后得到一个满意的穿透值。 

FKN的收敛性要求和穿透太大产生的计算误差总会是一对矛盾。解决此矛盾的办法是在接触算法中采用扩展拉格朗日乘子法。此方法在接触问题的求解控制中可以有更多更灵活的控制。可以更快的实现一个需要的穿透极限。 

2．拉格朗日乘子法与扩展拉格朗日乘子法 

拉格朗日乘子法与罚函数法不同，不是采用力与位移的关系来求接触力，而是把接触力作为一个独立自由度。因此这里不需要进行迭代，而是在方程里直接求出接触力（接触压力）来。 

Kx=F+Fcontact 

从而，拉格朗日乘子法不需要定义人为的接触刚度去满足接触面间不可穿透的条件，可以直接实现穿透为零的真实接触条件，这是罚函数法所不可能实现的。使用拉格朗日乘子法有下列注意事项： 

1）刚度矩阵中将有零对角元，使有些求解器不克使用。只能使用直接法求解器，例如波前法或系数求解器。而PCG之类迭代求解器是不能用于有零主元问题的。 

2）由于增加了额外的自由度，刚度阵变大了。 

3）一个可能发生的严重问题，就是在接触状态发生变化时，例如从接触到分离，从分离到接触，此时接触力有个突变，产生chattering（接触状态的振动式交替改变）。如何控制这种chattering，是纯粹拉格朗日法所难以解决的。 

因此，为控制chattering，ANSYS采用的是罚函数法与拉格朗日法混合的扩展拉格朗日乘子法。在扩展拉格朗日法中，可以采用实常数TOLN来控制最大允许穿透值。还有最大允许拉力FTOL。这两个参数只对扩展拉格朗日乘子法有效。 

在扩展拉格朗日乘子法里，程序按照罚函数法开始，与纯粹拉格朗日法类似，用TOLN来控制最大允许穿透值。如果迭代中发现穿透大于允许的TOLN值，（对178单元是TOLN，而对面面接触单元171-174则是FTOLN）则将各个接触单元的接触刚度加上接触力乘以拉格朗日乘子的数值。因此，这种扩展拉格朗日法是不停更新接触刚度的罚函数法，这种更新不断重复，直到计算的穿透值小于允许值为止。 

尽管与拉格朗日法相比，扩展拉格朗日法的穿透并不是零，与罚函数法相比，可能迭带次数会更多。扩展拉格朗日法有下列优点： 

1）较少病态，个接触单元的接触刚度取值可能更合理。 

2）与罚函数法相比较少病态，与单纯的拉格朗日法相比，没有刚度阵零对角元。因此在选择求解器上没有限制，PCG等迭代求解器都可以应用。 

3）用户可以自由控制允许的穿透值TOLN。（如果输入了TOLN，而使用罚函数法，则程序忽略它）