Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)

abaqus用户子程序fric,在接触分析中,定义复杂的摩擦模型,或者在热力耦合分析中,定义摩擦生热时,潜力巨大。这里先将子程序相关的基础知识,进行了整理。后续会更新基于子程序的相关应用案例。

1.概述:

用户子程序FRIC对应于关键字*FRICTION(定义一个摩擦模型。用于将摩擦特性引入表面接触模型中,来控制接触表面、接触对或连接器单元的切向接触行为。),以及交互界面里的接触属性中切向行为的所有内容(除了用户自定义外,abaqus中可以定义5种类型的摩擦行为(摩擦公式),每个公式中,主要是定义三方面的内容:摩擦因子,剪切应力,弹性滑动(可以恢复的滑动位移))。

Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图1Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图2Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图3Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图4

Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图5Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图6

用户子程序FRIC:

  • 可用于定义接触面间的摩擦行为;

  • 当Abaqus中提供的经典库仑摩擦模型的扩展版本限制太严格,或者需要在接触面间定义更复杂的切向应力时,可以使用;

  • 当接触属性模型包含用户子程序定义的摩擦时,当接触点闭合时,接触对的从属表面上的节点或者接触单元的积分点会调用子程序;

    每个增量步里的每次迭代,接触对中,从表面上,处于接触闭合状态的节点,会调用子程序。

  • 必须提供接触切向行为的完整定义;

  • 可以用来更新解相关的状态变量(例如常见的:弹性能量密度,摩擦生热等);

2.子程序基础介绍

子程序接口

SUBROUTINE FRIC(LM,TAU,DDTDDG,DDTDDP,DSLIP,SED,SFD,
1 DDTDDT,PNEWDT,STATEV,DGAM,TAULM,PRESS,DPRESS,DDPDDH,SLIP,
2 KSTEP,KINC,TIME,DTIME,NOEL,CINAME,SLNAME,MSNAME,NPT,NODE,
3 NPATCH,COORDS,RCOORD,DROT,TEMP,PREDEF,NFDIR,MCRD,NPRED,
4 NSTATV,CHRLNGTH,PROPS,NPROPS)
C
INCLUDE 'ABA_PARAM.INC'
C
CHARACTER*80 CINAME,SLNAME,MSNAME
C
DIMENSION TAU(NFDIR),DDTDDG(NFDIR,NFDIR),DDTDDP(NFDIR),
1 DSLIP(NFDIR),DDTDDT(NFDIR,2),STATEV(*),DGAM(NFDIR),
2 TAULM(NFDIR),SLIP(NFDIR),TIME(2),COORDS(MCRD),
1.1.8–1
FRIC
3 RCOORD(MCRD),DROT(2,2),TEMP(2),PREDEF(2,*),PROPS(NPROPS)
user coding to define LM, TAU, DDTDDG, DDTDDP, and, optionally, DSLIP, SED, SFD, DDTDDT, PNEWDT, STATEV
RETURN
END

必须要定义变量:

1. LM

相对滑动标记。仅在确定接触点的接触状态处于closed时,用户子程序FRIC才会被调用;即两种情形:接触压力是正的,接触点在前一次迭代中,接触状态是closed;接触点处于过度闭合状态(干涉),接触点在前一次迭代中处于open状态。

在迭代过程中,LM的值会首先传递到子程序,此时的值为前一个迭代中定义的值。在一个增量开始时,或者如果接触点在前一个迭代中是open状态,这个变量将是否传递进子程序中,取决于前一个增量中的接触条件。如果接触点滑动,LM = 0;如果接触点是粘着的,LM = 1;如果接触点是open状态,LM等于2。

如果这次迭代中,相对滑动是允许的(由于滑移或弹性粘着),设置LM等于0。在这种情况下,子程序中必须将切向摩擦应力定义为相对滑动、界面接触压力p和其他预定义或用户定义的状态变量的函数。此外,子程序必须定义切向摩擦应力相对于接触压力p的倒数(即切向应力的变化与接触压力的变化之间的关系)。

如果这次迭代中,不允许相对滑动,设置LM等于1。界面处的刚性粘着条件,采用拉格朗日乘子方法施加。在这种情况下,不需要更新其他变量。如果LM总是设置为1,会生成一个“非常粗糙”的界面。当使用有限滑动,surface to surface的接触公式时,不建议将LM设为1。

如果摩擦被忽略(假设无摩擦的滑动),在这种情况下,不需要更新其他变量。如果LM总是设置为2,则生成一个“完美光滑”的界面。

也可以根据增量滑动量的信息,以及计算得到的摩擦应力来决定粘着/滑动状态(摩擦力计算得到摩擦应力,然后求得等效摩擦应力,去与临界剪切应力比较)。这些变量也会从abaqus/standard中传递到子程序中。

为了避免通常的摩擦接触问题中出现收敛问题,如果在前一个增量结束时,该接触点已打开,将LM设置为2,并且退出子程序,也就是说,如果Abaqus/Standard 设置 LM=2,当它调用这个子程序时,只需退出这个子程序。

Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图7

2. 如果LM返回的值是0

数组TAU(NFDIR)

在增量开始时,这个值传入子程序,其值是摩擦切向应力分量,必须在增量结束时更新这个值。

将等效剪切应力与临界剪切应力相比较,可以判断时粘着,还是滑动。

Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图8Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图9

二维数组DDTDDG(NFDIR,NFDIR)

方向b上的,摩擦切向应力相对于方向a上的相对运动的偏导数。

二维数组DDTDDP(NFDIR)

方向a上的摩擦应力对接触压力的偏导数。由于这些项对刚度矩阵产生非对称影响,只有在使用非对称方程求解器时才使用它们。

可以被更新的变量:

DSLIP(NFDIR)

不可恢复滑动位移(滑移)的增量。如果在先前的增量步中,LM是0,这个数组传递进子程序,其值作为前一次迭代中用户定义的值。否则,其值将是0。仅仅当LM的返回值是0时,这个数组才会被更新。

这个数组用于检测迭代之间的滑动反转。输出选项使用它来指示该点是粘着状态还是滑动状态。当一个增量步收敛时, DSLIP(NFDIR)中的值会在累积在SLIP(NFDIR)中,作为塑性应变存储。

Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图10

如图,一开始不发生滑动(粘着),随着摩擦剪切力增大,到达临界剪切应力(u*p),开始滑动,滑动后,随着滑动速度降为0,摩擦剪切应力会降低。从而由滑动又变成粘着状态。随着变形的进行,粘着后,还会发生反向滑动。

Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图11

其中t是滑动方向,这是由接触公式决定的(N-S,S-S,小滑动,有限滑动)。U是位移的增量。例如一个圆盘的旋转。U=速度(ω*r)*Dtime(增量步的时间).

Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图12

Abaqus子程序系列:FRIC(定义接触表面的摩擦行为)的图13

SED

该变量在增量开始时作为弹性能量密度传入,并应在增量结束时更新为弹性能量密度。此变量仅用于输出,对其他解变量没有影响。

SFD

该变量应定义为摩擦耗散增量。如果调用FRIC的接触单元或接触对称使用应力而不是力,单位是每单位面积的能量。对于常规应力分析,该变量仅用于输出,对其他解变量没有影响。在温度-位移耦合和热-电-结构耦合分析中,如果采用间隙生热模型,可将耗散转化为热量(也是fric子程序,比较重要的应用方向)。如果SFD没有定义,可以通过滑移增量DSLP、摩擦切向应力TAU的乘积得到的耗散来计算生热(生热(摩擦耗散、摩擦做功)=摩擦切向力(=摩擦因子*法向压力)*滑动位移)。

DDTDDT(NFDIR,2)

a方向上的摩擦切向应力,对两个表面温度的偏导数。这只对于温度-位移耦合和热电结构耦合单元是需要的,这些单元中,摩擦剪切应力是表面温度的函数。

PNEWDT

建议的新的时间增量步尺寸与现在使用的时间增量步尺寸的比值。如果使用自动时间增量步算法,会允许这个变量。

每次调用FRIC子程序之前,这个值会设置成一个非常大的值。

如果PNEWDT被重新定义,小于1,abaqus/standard将放弃此时的时间增量步尺寸,尝试用更小的时间增量尺寸。为自动时间积分算法提供的建议的新时间增量等于PNEWDT乘以DTIME(此时的时间增量步尺寸),其中使用的PNEWDT��所���调用的用户子程序(该用户子程序允许在此迭代中重新定义PNEWDT值)里定义的最小值。

如果对于这个迭代中,所有调用的子程序中,PNEWDT是一个给定的远大于1的值,并且这个增量步在这个迭代收敛,abaqus/standard会增加时间增量步尺寸。为自动时间积分算法提供的新的时间增量为PNEWDT *DTIME,其中使用的PNEWDT是所有调用的用户子程序(该用户子程序允许在此迭代中重新定义PNEWDT值)里定义的最小值。

如果在分析过程中没有选择自动时间增量算法,PNEWDT大于1.0的值将被忽略,小于1.0的PNEWDT值将导致作业终止。

STATEV(NSTATV)

包含用户定义的解相关状态变量的数组。可以指定可用状态变量的数量。这个数组将在增量步开始时传入,其中包含这些变量的值。如果将任何与解相关的状态变量与摩擦行为一起使用,则必须在增量结束时,在此子程序中更新它们的值。

仅仅为了传递信息进子程序的变量

DGAM(NFDIR)

如果在先前的迭代中,LM设置为0,它的值是这个增量步中的滑动增量。否则,为0。通过与DSLIP(NFDIR)进行比较,可以确定是否在此时,这个点由粘着变成滑动状态,或者在此时是否发生了滑动方向反转。

TAULM(NFDIR)

如果先前的迭代中,LM设置为1,该值是增量步结束时约束应力(等效剪切应力)的当前值。否则,将是0。通过与临界剪应力的比较,可以确定此时这个点是否由粘着变成了滑动。

有些量是对于增量步更新的,有些量是对于增量步中的每次迭代更新的。如上面的DGAM和TAULM是每次迭代中更新的,而DSLIP是增量步中更新的。

PRESS

增量步结束时的接触压力

DPRESS

接触压力的增量

DDPDDH

在软接触的情况下,此时的接触刚度

SLIP(NFDIR)

在增量步开始时,总的不可恢复滑动。这个值是从前面的增量步中,DSLIP(NFDIR)的累计值。

KSTEP

分析步号

KINC

增量步号

TIME(1)

增量步结束时,分析步的时间值

TIME(2)

增量步结束时,总时间的值

DTIME

这个增量步的时间尺寸

NOEL

接触单元的单元编号。如果定义的是接触表面,传递进来的是0

CINAME

与摩擦定义关联的用户定义的表面的名称,左对齐。对于接触单元,它是与摩擦定义相关联的界面定义所给定的单元集名称。如果给界面定义分配了一个可选名称,那么CINAME将作为这个名称传入,左对齐。

SLNAME

从表面名称。如果使用接触单元,传入空格。

NSNAME

主表面名称。如果使用接触单元,传入空格。

NPT

接触单元的积分点数量。如果接触表面被定义,传入0。

NODE

与这个接触点相关的用户定义的全局从节点号(或根据部件实例的组装定义的模型的内部节点编号)。如果使用面对面接触公式,则对应于约束的主要从节点。如果从接触单元调用,则作为零传递。

COORDS(MCRD)

包含这个点,此时坐标的数组

RCOORD(MCRD)

如果主表面被定义为刚性表面,则传递此数组,其中包含刚性表面上对应点的当前位置和方向的坐标。

DROT(2,2)

旋转增量矩阵。对于与三维刚性表面的接触,这个矩阵表示相对于刚性表面的表面方向的增量旋转。这样做是为了在这个子程序中,将向量值或张力值状态变量进行正确的旋转。在调用FRIC之前,应力和滑动分量已经旋转了这个量。该矩阵作为单位矩阵传入二维和轴对称接触问题。

TEMP(2)

从节点和对应的主表面,此时的温度

PREDEF(2,NPRED)

一个数组,包含当前增量步结束时的所有用户指定场变量的一对值(分析开始时的初始值和分析期间的当前值)。从接触对调用FRIC,一对值中的第一个值对应于从节点,第二个值对应于主表面上最近的点。如果从一个大滑动接触单元调用FRIC,则PREDEF(1,NPRED)对应于该单元积分点处的值,PFREDEF(2,NPRED)对应于对应表面上的最近点。如果从一个小滑动接触单元调用FRIC,则PREDEF(1,NPRED)对应于第一侧积分点处的值,PFREDEF(2,NPRED)对应于该单元对应面上积分点处的值。

NFDIR

摩擦方向的数量

MCRD

接触点处坐标方向的数量

NPRED

预定义场变量的数量

NSTATV

用户定义状态变量的数量

CHRLNGTH

接触表面的面特征尺寸,可以用于定义最大允许的弹性滑动

PROPS(NPROPS)

用户定义的属性值数组,用来定义接触表面间的摩擦行为

NPROPS

用户定义的,与摩擦模型相关的属性值的数量

(13条)
默认 最新
我使用FRIC子程序静态接触分析,结果应力和接触云图和abaqus内置模型对的上,但是模态分析结果却显示两个接触的物体振型没连上。请问FRIC子程序静态接触分析后能否进行模态分析呢?谢谢
评论 点赞
谢谢分享,学习了。
评论 点赞

查看更多评论 >

点赞 35 评论 12 收藏 34
关注