基于微观接触模型的支撑基座刚度分析

曹顺成,张志强,肖 博,孙韵韵

(武汉大学动力与机械学院,武汉 430072)

0 引言

在中国各地同步辐射光源中,电子储存环是主体,为了向用户提供在空间位置上高度稳定的同步辐射光,电子束轨道的稳定需要被控制在微米量级,且支架系统振动本征频率需要足够高,避免引起共振[1]。目前多数同步辐射光源支架系统采用一种楔形运动机构[2]作为支撑基座以完成位置升降功能,董玉玺[3]提出支撑基座不同分布方式对支架系统的本征频率影响较大,ZHANG等[4]提出支架系统良好的刚度和阻尼能显著提高X射线束的稳定性。其中支撑基座的刚度计算作为支架系统模态分析的重要依据,缺少相关计算理论支撑和详细计算流程。

工程中最常用的手段是使用有限元软件进行刚度分析,各零件间的接触方式常采用共节点或接触对[5],并把整体机构等效为一个弹簧,根据胡克定律计算整体刚度。这种宏观结构模型可以简化计算规模,但不够精确。在一些关键部件的实际接触中,压力相对较高,导致表面发生塑性变形,宏观结构模型有限元分析会忽略屈服效应造成的影响,导致刚度计算值大于实际值,使支架系统本征频率计算值偏高。

通过微观粗糙表面处理宏观接触表面,可建立更精准的动力学模型,WANG等[6]使用微观粗糙表面表达式建立螺栓机械接合面的动力学模型,和模态试验结果对比验证,具有较高的精度。为研究微观粗糙表面接触特性,已有大量学者建立了不同接触模型,比较经典的有GW模型[7]、ZMC模型[8]、CEB模型[9]和KE模型[10-11],理论接触模型受限于多种简化假设,计算过程复杂且不够直观,而有限元法可还原实际表面纹理信息,并充分利用计算机能力,处理复杂问题较灵活。

针对支撑基座的竖直刚度计算,现存问题是有限元宏观计算结果不够精准,表面间的接触状态不清晰,影响刚度的原因不明朗,本文提出一种等效刚度模型,其变形更贴近实际接触效果,研究支撑基座在两种分析方式下的刚度差异和变化规律,探究影响刚度的因素,有助于优化结构设计,对提升电子储存环支架系统稳定性有重大意义。

1 微观粗糙表面有限元建模方法

1.1 基本理论

本文研究的支撑基座如图1所示,其中上法兰与支架底面贴合,向下传导载荷;关节轴承作为万向结构,允许在调节时任意方向转动一定角度从而使承载面与目标底面紧密贴合;楔形块机构可将下楔块横向位移转为上楔块高度升降,布置在混凝土地基上。

图1 支撑基座实物工作图

广义地说,任何机械结合面可看作弹性系统,可利用有限元子结构法思想,将宏观接触表面划分为若干个微观粗糙表面,先分别确定各子结构刚度特性,再将子结构装配成整体结构,获得组合结构等效刚度,其动力学模型如图2所示,每个接触单元视为弹簧元件,则结合面的刚度模型可写为:

图2 机械结合面的动力学模型

(1)

式中:K是宏观接触表面刚度,kn是微观粗糙表面刚度,S/Sn是宏观接触表面和微观接触单元面积比值,该式称为刚度并联公式。计算出一个微观接触单元的刚度,即可推导出宏观接触面的刚度,因此对微观粗糙表面的分析需要比较精确,重构模型的还原度要求高。

1.2 重构随机粗糙表面

为了使重构模型更贴近实际轮廓,需要获取原始粗糙表面的微观纹理信息。本文数据采集使用的测量设备是NanoFocus非接触式激光扫描轮廓仪,材料为45#钢,通过铣削、车削和磨削工件表面,获得3种加工方式下的表面形貌,如二维纹理特征图3所示。

(a) 铣削 (b) 车削(c) 磨削

大量研究结果发现,粗糙表面上微凸峰的高度分布满足高斯分布规律,其表面轮廓具有指数自相关函数关系。在陈辉、吴少雷等[12-13]提出的粗糙表面模拟方法的基础上,本文利用时间序列模型、滤波技术模拟生成了具有指定自相关函数的粗糙表面,图4为粗糙表面的模拟步骤。

图4 随机粗糙表面重构步骤

根据上述随机粗糙表面的生成方法,利用MATLAB软件构造模型和调节模拟表面的形貌参数,用于拟合铣削工件表面的截面轮廓,两者起伏曲线的对比图如图5所示,高度分布的标准差相同情况下,粗糙峰密度偏差不超过10%,模拟粗糙表面已满足实际分析需求。

(a) 真实工件表面轮廓 (b) 模拟工件表面轮廓

拟合数据可通过ANSYS命令直接生成三维点云模型,再采用自下而上的实体建模方式(从点到线、从线到面、从面到体)生成图6中的有限元微观粗糙表面模型。

图6 有限元微观粗糙表面模型

2 微观尺度下的接触分析

2.1 粗糙表面接触有限元分析

根据MCCOOL[14]提出两个粗糙表面的接触可等效转换为一个光滑的刚性平板和另一个粗糙的弹性板块相接触。两个粗糙表面的材料原为45#钢,则转换后弹性板块的弹性模量E=2.0×105MPa,泊松比μ=0.26,屈服强度为σs=355 MPa,延伸率为δ=16%,平板刚度足够大保证不变形。

本文通过仿真计算弹性板块粗糙表面的应力应变关系,推导出接触刚度。基体采用实体单元(SOLID187单元)进行网格离散和模拟,随后创建接触对,刚性平面作为目标面(TARGE170),粗糙表面当作接触面(CONTA174)。接触面网格最密集,沿着厚度方向渐变式减少网格数量,最终节点数量约为56万,计算规模可观。

刚性平板侧面采用无摩擦支撑约束,仅能竖直方向运动,弹性板块底面采用固定约束,用刚性平板给粗糙表面施压,为了保证模型计算能够收敛,初始提供较小的位移确保刚性平板与粗糙平面之间接触生效,并利用大变形和弱弹簧功能加以辅助。

将宏观接触表面的压力之和转换成单个微观接触单元的接触力,需要根据面积比例缩减载荷。对面积为4 mm2的粗糙表面施加载荷,其受力变形结果如图7所示,随着法向载荷增大,产生变形的范围不断增大,这是因为部分微凸体高度下降,结合面的接触间隙减小,微凸体接触数目增多,粗糙表面间的真实接触面积随之增大。

图7 粗糙表面受力变形图

2.2 验证有限元接触模型计算

为验证ANSYS软件仿真微观粗糙表面的数据合理性,本文列举GW模型、ZMC模型、CEB模型、KE模型的理论曲线与ANSYS有限元仿真曲线进行比较。理论模型的数值计算必须代入参数,而粗糙峰的平均曲率难以获得,根据NURI等[15]进行经典工程表面试验测得的数据,结合实际轮廓,得到面密度η=6.209×10-3与仿真模型相近,选用参数如表1所示。

表1 仿真参数

利用2.1节计算的粗糙表面受力变形结果,根据胡克定律可计算出法向刚度。由于微观粗糙表面尺寸仅达毫米级,相较于宏观表面可忽略不计,于是将该区域刚度近似为某一点刚度。

分别按3种加工方式(铣、车、磨)的形貌构建模拟粗糙表面和理论计算参数,将不同接触模型下结合面法向接触刚度随载荷的变化曲线进行对比,如图8所示,可用于验证本文的微观三维模型接触刚度分析的有效性和可靠性。

(a) 铣削 (b) 车削

由图8可知,ANSYS模型的法向刚度与GW模型、ZMC模型、KE模型以及CEB模型的变化趋势是基本一致的,从3幅图中可以发现以下几点。

(1)随着法向载荷增大,接触刚度逐渐增大,超过一定载荷后,增长速率逐渐减小,因为部分微凸体进入塑性屈服阶段,导致局部接触刚度降低,且接触面积受粗糙峰高斯分布函数的影响,到达一定程度后趋于饱和,所以接触刚度变化趋于平缓。

(2)不同模型在相同载荷下的接触刚度有区别。其中,ANSYS模型与KE模型的曲线最接近,这是由于两者都采用有限元法研究接触问题,更贴近真实接触状况。相似的曲线变化规律进一步说明了ANSYS接触刚度模型的有效性。

3 有限元宏观模型和微观表面刚度差异

3.1 支撑基座等效刚度

通过微观接触分析可获得微观粗糙表面的刚度变化曲线,而本文研究对象是支撑基座,如图9所示,多个接触面组成了整个支撑基座的弹性系统,由微观表面刚度计算支撑基座竖直刚度的步骤如图10所示。先根据宏观模型有限元分析结果,获得接触压力分布图,再利用刚度并联式(1)计算出每个宏观表面的接触刚度ki(i=A～E);随后由串联式(2)推导出支撑基座整体竖直刚度K。

图9 支撑基座结构剖面图

图10 支撑基座刚度计算流程

(2)

式(2)采用弹簧串联原理,刚度为ki(i=A～E)的接触面在同一轴线上串联,在外力的作用下所有结合面都被压缩,串联等效刚度为K。本文重点分析运动结合面(A、D)对刚度的影响,固定结合面(B、C、E)由于零件间固定联结,可看作一体,其刚度由材料弹性模量决定。

3.2 接触对的等效刚度计算方法

为了使用微观尺度刚度来建立具有宏观尺度的机械结合面的动力学模型,必须通过模拟来获得接触压力分布和实际接触面积。建立有限元宏观模型如图11a所示,约束各零件自由度,建立上楔块和下楔块接触对(图11b),建立关节轴承球碗和球头接触对(图11c),上法兰顶面建立刚性区域传递支架的压力,后处理中查看应力分布图可得知结合面接触区域面积。工程中支撑基座均以铸铁支架和重型磁铁[16]为承重目标,本节以负载为10 000 N为例,对两个接触对进行ANSYS仿真计算,采用罚函数算法处理接触状态,法向刚度比例因子FKN=1,算法将根据材料参数获得接触刚度估计值。

(a) 整体结构模型 (b) 楔形块接触对 (c) 关节轴承接触对

关节轴承采用GE 50ES型号,材料为轴承钢。由于其接触面是弧面,承载不均匀会导致表面变形,挤压后实际接触区域缩小,如图12所示,实际接触上边界与水平线夹角为θ1,下边界与水平线夹角为θ2,可见上边界存在偏移现象(理想的接触上边界应与水平线齐平),关节球头整体下沉。

图12 接触面D受力变形图

对曲面接触区域的竖直刚度分析采用积分累加方法,假设区域内各点的法向载荷相等,则法向接触刚度相等,刚度分解方式如图13所示,若弧面接触区域内某一点的法向刚度为kd,则其刚度垂直分量kd-n计算公式如式(3)。

图13 接触面D刚度分解图

kd-n=kd×sinθ

(3)

根据等效刚度并联式(1),可得到实际接触的环形区域在竖直方向上的有效刚度kD表达式:

(4)

式中:R是关节轴承球头的球半径,kD是通过各点的竖直刚度kd积分而来。

3.3 结果对比和分析

利用有限元法分析支撑基座竖直刚度的计算模型可分为两种:装配体仿真常用的宏观接触模型和本文提出的等效刚度模型,计算方法已在图10 给出。通过两种计算模型得到接触面A、接触面D和支撑基座整体结构的不同刚度曲线,研究两者变化规律的异同。

接触面A的竖直刚度kA是楔形面接触刚度的垂直分量,其曲线如图14所示,对表面施加0～16 000 N,两者的竖直刚度均随着载荷逐渐增大,刚度增长速率逐渐减小,宏观接触模型数值更大,其原因是粗糙表面考虑微凸体高度分布,表面间不能完美地贴合,而宏观模型的表面接触平滑,真实接触面积较大,导致刚度计算值大于实际值。载荷增大过程中,由于粗糙表面刚度考虑了塑形变形的影响,局部接触刚度降低,这种趋势也符合2.2节分析的结果。

图14 接触面A刚度对比图 图15 接触面D刚度对比图

接触面D的竖直刚度kD曲线如图15所示,两者增长速率有明显区别,其中有两点原因:其一,接触面D实际接触区域较小,接触压力较高,此时粗糙表面产生塑性变形的区域占比较高,对刚度的削减更严重;其二,宏观接触模型涉及基体变形,法向载荷通过表面传递到基体上,球碗基体高度下降,导致曲面弧度增大,更有利于竖直方向的支撑,而等效刚度模型是以初始弧度计算,刚度分解后削弱严重,相同面积下积分结果差别很大。

支撑基座整体结构的竖直刚度K曲线如图16所示,整体结构竖直刚度和单个接触面刚度的变化趋势相似,随着法向载荷的增大,刚度逐渐增大,刚度增长率逐渐减小,这是因为整体结构是由各个接触面串联组成,而整体刚度却比任何一个接触面的刚度都小。其中接触面D对整体结构的影响最大,这是因为接触面D刚度最小,整体刚度的大小和数量级受限于接触面D,两种计算模型下的结果差异也主要由接触面D造成。

图16 支撑基座刚度对比图

4 总结

通过数值模拟方法构建出具有自相关函数的微观粗糙表面,与轮廓扫描仪探测的工件真实表面轮廓进行拟合,将表面数据导入ANSYS软件中构建粗糙表面实体模型,并将有限元仿真结构与不同理论接触模型对比,验证其有效性。随后,基于子结构法提出等效刚度模型,通过分析支撑基座接触面的受力状态,得到真实接触区域,结合微观粗糙表面的刚度变化曲线,利用串并联刚度理论计算出组合结构的竖直刚度,与宏观结构模型刚度对比,分析产生差异的原因。主要结论如下:

(1)ANSYS接触模型刚度曲线与GW模型、ZMC模型、KE模型以及CEB模型的变化趋势基本一致,随着法向载荷增大,接触刚度逐渐增大,超过一定载荷后,增长速率逐渐减小,相似的曲线变化规律说明了ANSYS接触模型的有效性。

(2)宏观接触模型和等效刚度模型皆受关节轴承接触面的影响最大,相同静态负载下,等效刚度模型刚度竖直刚度更小,且负载越大,两者差异越明显,主要由接触面的塑形变形和曲面的基体变形所引起。