吕 昕
(南京航空航天大学机械结构力学及控制国家重点实验室 江苏南京 210016)
随着旋转机械在航空航天、低温制冷、石油化工和农业生产等领域中被越来越广泛地应用,人们对其中的关键支承部件——轴承也提出了越来越高的要求。传统滚动轴承和液体润滑轴承都因其固有局限性而难以向超高转速和超高精度发展,由此气体润滑轴承便应运而生。
气体润滑轴承因其润滑介质黏度受温度影响较小[1],不仅可以满足高转速、高精度的需求,还有耐寒耐热、低摩擦、无污染、使用寿命长等优点。但气体润滑轴承同时也存在启停磨损严重、易失稳、支承刚度小和抗外界冲击能力弱等问题。而引入弹性箔片恰好能解决上述问题。箔片动压轴承是以箔片结构作为弹性支承的自适应轴承,箔片结构可随受载调整气膜厚度以适应间隙内的压力变化。可见箔片动压轴承对制造精度和转子对中性要求降低,适应环境能力强[2-3],且自身阻尼及库仑摩擦力还能抵消轴承-转子系统的涡动能量并抑制自激振荡的产生。
国外对动压轴承的理论分析起步较早,研究范围也更为全面。1982年,HESHMAT等[4]对单叶片和多叶片弹性箔片轴承进行了对比实验研究,发现多叶片轴承各叶片均可形成动压收敛间隙从而具有变刚度特性。次年,HESHMAT等[5]将箔片结构等效为平/波箔整体位移的线性弹簧模型从而推导出Heshmat公式,填补了箔片结构等效模型的空白。1995年,HESHMAT和HESHMAT[6]初步研究了悬臂型动压轴承,假设初始气膜厚度为抛物线分布、弹性支承构件为弹簧模型,此举简化了气弹耦合计算过程。2008年,DELLACORTE等[7]探究了第一代和第二代波箔动压轴承静态特性的影响机制,总结归纳了两代轴承的设计方法、制造工艺以及试验研究。2012年,莫霍克创新科技公司(Mohawk Innovative Technology,Inc.)结合悬臂型轴承和波箔型轴承的特点,提出了新的多叶式波箔动压轴承:在相邻平箔片搭接区域采用了弹性箔片支承结构。2014年,GAD和KANEKO[8]建立了新型箔片轴承的刚度模型,采用卡氏第二定理对箔片变形和位移进行求解计算。同年,DU等[9]对于含支承波箔的多叶式动压轴承进行了多种模型等效研究,发现与Heshmat模型吻合程度较好。
20世纪90年代,弹性箔片轴承技术及其理论才被引进入中国,由于起步较晚加上国外技术封锁,导致我国的理论研究和工程应用都处于相对滞后阶段。HOU等[10]对刚性轴承和波箔型轴承进行了试验研究对比,发现波箔型轴承的稳定性相较刚性轴承更好。徐润和马希直[11]应用弹性薄壳单元模型,耦合薄膜应力和平板弯曲效应求解第一代波箔动压轴承的静态特性。2018年,皮骏等人[12]应用有限元法和松弛迭代法进行差分迭代求解交错式箔片轴承的静力学特性。但交错式箔片轴承仅改变波箔结构,而平箔仍为传统整体式。次年,阮琪等人[13]针对悬臂型动压轴承,分析了气固特性参数对其承载特性的影响。2020年,燕震雷和伍林[14]就可倾瓦动压轴承性能受稀薄效应的影响进行了研究。但他们所研究的动压轴承均并未引入支承波箔。同年,冯凯等人[15]对新型三瓣式径向箔片动压轴承进行了热特性分析研究,但这一新型轴承仍采用整体式平箔。
综上所述,自多叶式波箔动压轴承提出以来,国外学者已对其进行了详细的理论分析和试验研究,而国内无论是在理论、试验还是工程应用上都对其研究甚少。鉴于未来轴承高精度、高转速的发展趋势,多叶式波箔动压轴承必将拥有广阔的应用前景。因此,建立可靠理论模型、充分掌握其运行机制,研究结构简单、适合国内工业水平的多叶式波箔动压轴承具有重要的理论价值和工程意义。本文作者以多箔叠加弹性结构作为切入点,应用悬臂弯曲梁模型获得不同于前三代轴承的气膜厚度方程;将其耦合动压Reynolds方程应用有限差分法进行迭代求解,获得多叶式波箔动压轴承的静态特性。在此基础上,研究转速、偏心率以及工况温度的改变对其静态特性的影响,为之后的工程设计提供参考。
文中研究的多叶式波箔径向动压轴承属于第四代波箔动压轴承,它在结构上结合了悬臂型轴承和波箔型轴承的特点:顶层平箔两两相互搭接,波箔片为平箔片提供弹性支承。另外,这一结构轴承还拥有第四代波箔动压轴承的两大优势:大预紧力技术和耐高温、耐磨涂层技术的应用。这些特点使得轴承的支承能力可达到刚性表面下的支承强度,耐高温、耐磨的涂层技术也为轴承高速运转产生的气动热问题提供了解决方法。
图1所示为多叶式波箔动压轴承结构示意图。
图1 轴承结构
该径向动压轴承的平/波箔叠加弹性结构由5片顶层平箔和5片支承波箔组成,顶层平箔相互两两搭接在一起,支承波箔安装在顶层平箔下,平箔键块安装在轴承套的键槽内固定,波箔一端的翻边固定在平箔键块和轴承套键槽内。其中,顶层平箔由箔片和键块组成,箔片和键块通过点焊焊接在一起。
当转子高速旋转时,转轴和箔片结构间会形成楔形间隙,随之生成气膜提供向上的支承力来平衡轴承所受载荷,从而实现转轴的悬浮旋转。为求解气膜压力,需要先分析非均匀分布的气膜厚度。考虑到轴承内部气体流动情况,先作以下假设:忽略气膜在层与层间滑移作用;忽略气膜压力沿气膜厚度方向变化;忽略转子曲率半径对气膜方向和形状影响;将黏性气体沿转子接触面相对运动视为平移运动,其速度大小等同于转子接触面切向速度;忽略惯性力及体积力作用;假设气体黏度和密度沿气膜厚度方向不变化。
在上述假设的基础上,可以写出可压缩气体的Reynolds方程:
(1)
式中:x为轴承周向坐标;z为轴承轴向坐标;h为气膜厚度(m);p为气膜压力(Pa);μ为气体动力黏度(Pa·s);U为转子沿周向运动速度(m/s)。
当动压轴承的散热特性良好时,可将气体动力黏度视为一个常量,也无需考虑温升对润滑气膜压力产生的影响。所以对定常气体,可将其Reynolds方程进行量纲一化:
(2)
其中量纲一化参数如下
(3)
式中:C为间隙厚度(m);L为轴承厚度(m);R为转子半径(m);ω为转子角速度(rad/s);pa为环境压力(Pa);θ为量纲一化周向坐标;λ为量纲一化轴向坐标;H为量纲一化气膜厚度;P为量纲一化气膜压力;Λ为轴承数,用以反映轴承的运行条件和性能参数。
为求解量纲一化的Reynolds方程,需要对气膜厚度方程进行推导。但在推导公式之前,考虑到箔片结构的实际情况,假设箔片结构刚度分布均匀并忽略平箔变形,只考虑其随波箔的位移。
由于支承波箔被等效为线性弹簧支承,不考虑顶层平箔刚度特性和内部摩擦阻尼作用,对箔片结构进行受力分析,可以写出量纲一气膜压力与量纲一化箔片结构变形量间的作用关系式:
(4)
式中:s为波箔单元长度(m);kb为波箔等效刚度(N/m);u为箔片结构变形量。
对平/波箔耦合的多箔叠加弹性复合结构进行受力分析时,由于5个相同结构的顶层平箔在轴承套内有序排列,故取其中之一进行分析即可。将平箔片视为悬臂弯曲梁模型而波箔片视为弹性支承模型。图2所示为平箔片等效的悬臂弯曲梁模型。图中点O为转子圆心,点O1为平箔片曲率圆心,点O′为轴心,点A为平箔片固定端,点B为平箔片自由端,点Cp为两相邻平箔片搭接点,点Tp为平箔片切圆切点,点p为平箔片上任意点。
图2 等效悬臂弯曲梁模型
xOy坐标系与x1Oy1坐标系之间的转换角度为
(5)
而如图3所示,任意坐标系xiOyi均可由动坐标系xOy旋转一个角度γ变换得到,T(γ)就是动坐标系xOy变换到任意坐标系xiOyi的变换矩阵:
图3 动坐标系与任意坐标系转换关系
(6)
令γ=-Ψ,在已知偏心距en和偏位角φn的情况下即可推得xOy坐标系下平箔片圆心O1、轴心O′和平箔片上任意点p的坐标:
(7)
(8)
(9)
令γ=Θ,i=1、2、3、4、5为各平箔片的编号,则各平箔片的圆心Oi和各平箔片上任意点pi的坐标可表示为如下形式:
(10)
(11)
计算气膜厚度的向量ai和bi的表达式如下:
(12)
即可求得平箔片上任意点pi处的气膜厚度:
(13)
为得到气膜厚度及压力分布情况,需要对可压缩气体的定常Reynolds方程进行求解,从而能进一步获得承载力和偏位角。在之前的研究中常用无限宽或短轴承理论来近似求解Reynolds方程,但这一方法精度较低,难以满足越来越精确的数值计算要求。所以,可在较短时间内获得良好数值计算结果的有限差分法现如今被大量应用于Reynolds方程的数值分析求解中。
有限差分法究其根本是导数的差商表示法,导数的差商表示法又分为前差商法、后差商法和中心差商法3种,这3种差商表达形式以中心差商法的计算精度最高,所以常被用于对Reynolds方程的离散化。其表达式如下所示:
(14)
式中:yj+1/2为点yj前半点处的值;yj-1/2为点yj后半点处的值;δ为计算步长。
在进行Reynolds方程离散化之前,需要对轴承内气膜进行网格划分,用各个节点处的压力值作为各阶差商的数据,从而将Reynolds方程离散成代数方程组。通过求解各个节点处的压力值,即可近似求得气膜压力分布情况。气膜网格划分如图4所示,沿θ方向均匀划分了m格,i编号从1到m+1,步长为Δθ=2π/m;沿λ方向均匀划分了n格,j编号从1到n+1,步长为Δλ=2π/n。
图4 气膜网格划分
由此可对量纲一化Reynolds方程中各项进行离散化后重新代入得到用节点(i,j)周围4个节点压力值来计算中央节点压力值的表达式:
(15)
Ei,j=ΔθΛ(Pi+1/2,jHi+1/2,j-Pi-1/2,jHi-1/2,j);
Fi,j=Ai,j+Bi,j+Ci,j+Di,j。
图5 耦合迭代求解流程
气膜厚度和压力方程分别不断循环更新,最终迭代终止的收敛条件如下:
(16)
当耦合迭代求解结果收敛,即可得到稳态下轴承求解域中所有节点的气膜厚度和压力分布值。承载力W可分为沿偏心方向的切向力Wt和垂直于偏心方向的法向力Wn,将求解得到的气膜压力分别按偏心方向以及垂直于偏心方向进行积分,它们的量纲一化表达式如下:
(17)
由此可以求解得到偏位角Φ的数学表达式:
(18)
文中基于有限差分法对气膜厚度方程和可压缩气体Reynolds方程进行流固耦合迭代求解。为验证方法及程序的正确性,以文献[16]中刚性、波箔动压轴承为算例进行结果验证。图6(a)和6(b)所示分别为文献[16]和计算程序所得气膜压力分布曲线。对比两图所得求解结果一致,验证了所用方法和自编程序的正确性。
图6 计算程序验证
量纲一承载能力、气膜压力及厚度分布、偏心率以及偏位角这些轴承参数都能用来表征动压轴承稳定运转时的静态性能。以多叶式波箔动压轴承为研究对象(其具体参数见表1),同时考虑气体的压膜效应与支承波箔的变形效应,采用有限差分法耦合求解Reynolds方程和量纲一气膜厚度方程,从而获得轴承的静态特性,进而总结多种气固特性参数对多箔径向动压轴承气膜压力、承载力以及偏位角的影响规律,为多叶式波箔径向动压轴承的模型设计、应用条件提供理论参考。
表1 轴承具体参数
在偏心率为0.5、转速为9×105r/min的工况条件下,可以得到多叶式波箔动压轴承的全尺寸气膜压力分布,如图7所示。可以很明显看出有5个压力峰,这分别对应了5片分离式平箔上的最大气膜压力。另外由于转子偏心不对中,导致其中一片顶层平箔上的最大气膜压力高于其余平箔片。
图7 轴承全尺寸气膜压力分布
图8所示是偏心率为0.7、不同转速条件下轴承中截面的气膜压力分布。可知,随着转速的升高,各箔片上的气膜压力都获得了不同程度的增大,其中受转子不对中影响最大的箔片上气膜压力提升最多。这是由于随着转子转速的提升,楔型间隙内黏滞气体会随转子一同高速旋转,密度减小,不可压缩性增大,从而能提供更大气膜压力。图9所示是转速为7×105r/min、不同偏心率条件下轴承中截面的气膜压力分布。可知,随着偏心率的增大,各箔片上的气膜压力同样都获得了不同程度的增大,其中也是受转子不对中影响最大的箔片压力提升最多。这是因为随着偏心率的增大,转子不对中性增强向某一方向偏移,从而导致该箔片上的气膜压力获得大幅提升而其余箔片上气膜压力趋近于1。
图8 不同转速下气膜压力分布
图9 不同偏心率下气膜压力分布
图10和11分别示出了承载力随偏心率和转速变化曲线。可知,承载力会随着偏心率或转速的增大而增大;其中承载力随转速增大呈线性增大而随着偏心率的增大呈抛物线式快速增大。究其原因是由于偏心率或转速的增大会导致各箔片上气膜压力不同程度的提升,相对应地作为反作用力,动压轴承承载力也会出现不同程度的增大。另外,对比转速、偏心率变化对气膜压力的影响可知,最大气膜压力、承载力变化受偏心率的影响要大于受转速的影响。
图10 承载力随偏心率变化曲线
图11 承载力随转速变化曲线
图12是不同转速下轴承偏位角随偏心率变化曲线。可知,偏位角会随着偏心率增大而快速减小,这与第一代波箔动压轴承的偏位角变化规律相一致。偏位角减小意味着轴承-转子系统偏移量减小,即随着偏心率的增大,轴承运转稳定性提升,这对轴承受扰后恢复稳定运转是有利的。
图12 偏位角随偏心率变化曲线
由于动压轴承的运转工况主要为高速轻载,所以常常伴随着大量的气动热。文中探究不同工况温度下轴承静态特性的变化。图13所示是不同工况温度下轴承气膜压力分布变化曲线。可知,随着温度升高,5片平箔上的气膜压力都获得了不同程度的提升。这是由于温升导致润滑气体的动力黏度增大,即润滑气体黏滞性增强,更易在楔型间隙内形成动压气膜,从而能提供更大的气膜压力。图14中最大气膜压力随温度变化曲线也印证了这一点:随着温度的升高,各个箔片上的最大气膜压力也都呈线性增大。因此,对于高速运转的箔片动压轴承而言,气动热问题是一个亟需解决的难题,如何做好轴承的冷却对于多叶式波箔气体动压轴承未来的应用具有重要意义。
图13 不同温度下气膜压力分布曲线
图14 最大气膜压力随温度变化曲线
(1)应用悬臂弯曲梁模型建立多箔叠加弹性结构的理论模型,基于有限差分法耦合气膜厚度方程和动压Reynolds方程进行迭代求解,得到多叶式波箔动压轴承的静态特性结果。
(2)研究了偏心率、转速的改变对多叶式波箔动压轴承静态特性的影响。结果表明,随着转速或偏心率的升高,各箔片上的气膜压力都获得了不同程度的增大,其中受转子不对中影响最大的箔片上气膜压力提升最多;承载力会随着转速增大呈线性增大而随偏心率增大呈抛物线式增大。
(3)初步探究了工况温度对轴承静态特性的影响,发现随着温度的升高,各个箔片上的最大气膜压力也都呈线性增大。因此,对于高速运转的多叶式波箔动压轴承而言,气动冷却是未来研究的热点与难点。