王相平,李 星,王剑锋,吴少培,丁旺才,李国芳
(1.兰州交通大学机电工程学院,甘肃 兰州 730070;2.包头铁道职业技术学院铁道机车车辆系,内蒙古 包头 014060)
中低速磁浮交通是中国具有自主知识产权的新技术,近年来得到了广泛发展。空气弹簧悬挂系统作为中低速磁浮车辆的关键隔振部件,其动力学性能直接影响车辆的安全性、舒适性,是决定中低速磁浮车辆发展前景的重要因素[1-2]。
近年来,以线性弹簧、阻尼并联的空气弹簧等效模型因其简单、高效的计算优势而被广泛应用于中低速磁浮车辆动力学模型中[3-6],但难以准确表征空气弹簧系统的幅频依赖性也致使其计算结果存在较大误差[7]。罗英昆等[8]运用AMESim 建立了考虑气动特性的空气弹簧非线性模型,并结合SIMPACK模拟了小半径竖曲线上高速磁浮车辆空气弹簧系统的动力响应;Berg[9]根据铁道车辆空气弹簧力学特征建立了包含弹性力、摩擦力和黏性力的空气弹簧一维非线性模型,该模型对求解精度和计算效率的合理折衷使得其在轨道车辆动力学分析中被广泛使用;戚壮等[10-11]进一步考虑空气弹簧胶囊材料的黏弹性特性,运用分数阶微积分理论修正、优化了铁道车辆空气弹簧气动热力学模型,并将其与UM 结合进行了整车动力学计算;陈俊杰等[12-13]采用摩擦模型和分数导数Kelvin-Voigt 模型对空气弹簧橡胶气囊的频率、振幅相关性进行描述,建立了汽车空气弹簧非线性模型。目前,国内外对中低速磁浮车辆空气弹簧动力学模型的相关研究极少,采用等效线性模型或将上述车辆空气弹簧动力学模型移植于磁浮车辆进行计算是最为常见的做法[1-8]。但中低速磁浮车辆无导向电磁铁,横向阻尼主要由空气弹簧悬挂装置提供[1],车辆横向欠阻尼特性导致其空气弹簧结构与一般铁道车辆存在显著差别。以中低速磁浮车用SRI1R160A 和SRI1R160B 系列空气弹簧[14-15]为例,与铁道车辆空气弹簧[16]相比,为满足其特定安装方式与功能特点,中低速磁浮车辆空气弹簧胶囊呈筒形并逐步演化出中央导柱。两者结构与功能上的巨大差异对现有铁道车辆空气弹簧模型是否适用于中低速磁浮车辆提出了新的质疑与挑战,在充分考虑二者结构特征的基础上,建立合理的中低速磁浮车辆空气弹簧动力学模型就显得极为必要。
鉴于此,首先在充分考虑中低速磁浮车辆空气弹簧结构特征的基础上,分析各向载荷作用下空簧的受力特征,建立各关键部件的力学模型,再依据力学等效原理,建立方便用于动力学计算的空气弹簧等效模型;其次进行系统参数辨识,并以某型空气弹簧为例对模型进行验证;最后将本文模型应用于中低速磁浮车辆,通过线路动态测试结果进一步验证本文模型工程应用的可靠性。研究结论以期为中低速磁浮车辆动力性能预测提供理论基础和技术保障。
全文总体技术路线如图1 所示。图1 中,F=表示以合力F等于各分力Fi之和建立空簧等效模型;线性刚度K由各部件刚度k1和k2叠加而成;阻尼C由节流孔阻尼效应等效计算得到,表示为与阻尼孔有效面积s相关的函数f(s)[8]。
图1 技术路线Fig.1 Technical route
中低速磁浮车辆空气弹簧悬挂系统由空气弹簧、附加气室、滑台、防脱开装置等部件组成,其安装位置及结构如图2 所示。图2 中,空气弹簧下部经锥形销定位后安装于悬浮托臂上的空气弹簧安装筒内,上部通过圆柱销定位后采用螺栓将其与滑台固接,滑台可在导轨上移动,导轨固接于车体,同时连接车体的还有位于悬浮架纵梁上的牵引拉杆。
图2 空气弹簧安装位置及结构示意图Fig.2 Installation position and structure diagram of air spring
空气弹簧系统结构如图3(a)所示。图3(a)中,H0为空簧标准工作高度;H1为上盖螺栓高度;H3为上盖螺栓距摩擦块距离;D0和D2分别为止挡和胶囊直径。由于空簧安装位置较小故将活塞制作为中空结构,使其同时兼具附加气室和应急止挡的作用。
图3 空气弹簧结构特征及其受力分析Fig.3 Structural characteristics and force analysis of air spring
1.2.1 垂向动力学模型
由于上盖板到摩擦块的距离H2小于下盖板到中央导柱底部的距离H4,因此垂向载荷Fz作用下空气弹簧垂向位移zi<H2时所受到的力主要有:胶囊本体变形产生的黏弹性力Fze、气动作用力Fza以及中央导柱与两侧壁面的摩擦力Fzf;当zi≥H2时上盖板与摩擦块接触,压迫应急止挡,故还需考虑应急止挡作用力Fzs,受力分析如图3(b)所示。垂向作用力表示为:
式中zu和zd分别为上盖板和活塞垂向位移。
胶囊由帘线橡胶复合材料制成,其变形产生的弹性力Fze为典型 的黏弹性力[10,12],故以分 数导数Kelvin-Voigt 模型表征;摩擦力Fzf与空簧垂向位移zi密切相关,故以库仑摩擦模型表征[9];应急止挡由橡胶与金属材料粘合而成,具有显著的非线性特性,故以非线性弹簧-阻尼模型等效;气动作用力Fza为空气弹簧系统主要的作用力,结合其几何特征依据气动热力学相关理论计算。各模型示意图如图4所示。
图4 各非线性模型示意图Fig.4 Schematic diagram of each nonlinear model
以分数导数Kelvin-Voigt 模型表征的胶囊黏弹性力Fze为[17-18]:
式中Kze为垂向线弹性刚度;bz为分数导数阻尼参数;az为分数导数阶次;Dazzi(t)表示位移zi(t)的az阶分数导数。
以库仑摩擦模型表征的摩擦力Fzf为[9]:
式中Ffs_z和zs分别为参考状态的初始力和初始位移;Ffmax_z为最大摩擦力;z2为时所对应的位移大小
以非线性弹簧-阻尼等效的应急止挡力Fzs为:
式中FCzs()和FKzs(zi)分别为阻尼力和弹性力;Kzs和Czs分别为其非线性刚度和阻尼。
气动作用力Fza为[13]:
式中Se(zi)为空气弹簧有效面积;p为胶囊内气体压强;pa为大气压强。
考虑中央导柱几何特征,膜式空气弹簧的有效面积Se(zi)为[19-20]:
式中κ为空簧高度对有效面积的影响因子[19-20];z0为空簧初始工作高度;Ra为压强变化后的气囊半径;RD为气囊初始半径;Tz为气囊厚度;δ=0.8;Q1和Q2为与铺层角相关的值;D1为中央导柱直径。
理想气体过程方程与体积表达式如下式所示:
式中p0为初始状态气体压强;V和V0分别为当前和初始状态气体体积;α为有效体积变化率;n为气体多变指数,取n=1.30~1.38。
结合式(10),可将式(6)中的气体压强p表示为下式形式:
将式(7)和(11)代入式(6),则气动作用力Fza为:
1.2.2 纵/横向动力学模型
空簧为对称结构,其纵/横向动力学模型一致[9]。本节纵/横向力求解的目的是构建其与垂向的力学/几何关系,如此仅完全解析某一方向作用力便可得到整体模型的力学关系。限于文章篇幅,后文仅给出横向模型计算过程,纵向计算与横向计算完全一致,仅将y替换为x即可,故不再赘述。
横向力Fy作用下空气弹簧上下盖板发生相对横移,此时系统所受作用力主要有:中央导柱横向力Fyf、胶囊体积变化引起的横向气动作用力Fya、胶囊自身变形产生的横向力Fye和应急止挡横向力Fys,受力分析如图3(c)所示。图3(c)中,yd和yu分别为活塞和上盖板横移量;yi为任意状态下的空簧横向位移。横向力Fy表示为:
应急止挡以金属材料为主,为各向同性材料,其横向力Fys可按式(5)计算。以空簧中轴线为Z轴,上盖板处为XOY平面建立坐标系,胶囊结构各向应力如图5 所示。
图5 胶囊结构应力分析Fig.5 Structural stress analysis of rubber balloon
在上盖板与应急止挡约束作用下,胶囊轴向变形横截面为平面上的环状结构,故轴向应力σφ为:
式中R0(zi),R1(zi)分别为与高度相关的胶囊内、外径。
在胶囊囊壁取径向长度dr、轴向长度单位1、周向dθ角度的微元体,如图5(c)所示。则由法线方向上∑F⊥=0 可得:
式中σr为径向应力;σθ为周向应力;r为所取微元体到回转中心的径向距离,r∈[R0,R1]。
周向应力σθ作用下囊壁周向纤维伸长,横截面圆环周长增大,半径对应增大,记增大量为w;径向应力σr作用下囊壁径向拉伸,胶囊壁厚增大,记增大量为dw,如图5(c)所示,则周向应变εθ、径向应变εr分别为:
结合式(17),对周向应变εθ求导,整理后可得:
根据应力-应变关系,εθ和εr可分别表示为:
式中E和μν分别为胶囊材料的弹性模量和泊松比。
联立式(18)和(19),可得胶囊变形的应力微分方程为:
系统的边界条件为:
结合式(14),由式(20)和(21)可得三向应力之间的函数关系式为:
则各向作用力为:
式中Si为i向应力作用面积,i分别表示轴向φ、周向θ和径向r。
则胶囊自身变形产生的横向力Fye为:
即最终可建立如式(24)和(25)所示的Fye与Fze(Fφ)的函数关系,在式(3)计算出轴向力Fze的基础上求得横向力Fye。
与垂向类似,横向气动作用力Fya为:
式中Se_y(yi)为影响横向刚度的有效面积,是空簧有效半径Re和气囊自由圆弧中心垂向位移Δzi的函数,可按下式计算[21]:
根据几何关系,圆弧中心垂向位移Δzi为[21]:
式中Ha为胶囊高度;Hb为气囊自由圆弧中心距胶囊横向边界的距离。Ha,Hb和H3等参数见图3(a)。
中央导柱横向力Fyf主要为导柱抵抗变形的作用力,其值较大,故以弹性模量表征刚度。
1.2.3 模型等效
依据各作用力特性所建立的空气弹簧垂向、横向动力学模型如图6(a),(b)所示。模型复杂、参数较多导致其难以直接用于车辆动力学计算,故基于力学等效原理,建立空气弹簧等效模型如图6(c)所示。
图6 空气弹簧动力学模型Fig.6 Dynamics model of air spring
图6(c)中,Kz和Cz分别为垂向刚度和阻尼,其各项组成如图6(a)所示;Ky和Cy分别为横向刚度和阻尼,其各项组成如图6(b)所示。图6(a),(b)中,符号K和C分别表示刚度和阻尼,其下标首字母z和y分别表示垂向和横向,下标第二位字母e,f,a和s分别表示弹性项、摩擦项、气动项和应急止挡非线性项。各向刚度按下式计算:
系统中的库仑阻尼、结构阻尼以及其他非黏性阻尼,均采用能量法转化为等效黏性阻尼,即:
式中 ΔE为黏性阻尼在一个周期内消耗的能量,库仑阻尼ΔE=-4μFzA,结构阻尼ΔE=-vA2,其中μ和v分别为摩擦系数和比例系数;A和ω0分别为激振力的振幅和圆频率。
需要辨识的系统参数有:库仑摩擦模型中的初始力Ffs_z和初始位移zs,最大摩擦力Ffmax_z及对应的位移z(2式(4));分数导数Kelvin-Voigt 模型中的刚度Kze、阻尼参数bz和导数阶次az(式(3));应急止挡、胶囊材料、气动作用力所对应的刚度与阻尼(图6)。库仑摩擦模型参数通过空气弹簧静态试验辨识;分数导数Kelvin-Voigt 模型参数通过动态试验辨识;刚度与阻尼通过滞回曲线辨识。力-位移滞回曲线如图7 所示[9]。
图7 力-位移滞回曲线Fig.7 Hysteretic curve of force-displacement
主气室静态试验加载频率极低,胶囊黏弹性特性可以忽略,此时主气室静刚度Kz仅取决于摩擦模型和气动特性,即[9]
其中,Kze+Kza由滞回曲线zi=±z0两点处斜率求得;KAB由静态滞回曲线两端点所连直线AB的斜率求得。
则库仑摩擦模型中相关参数之间的关系为:
式中Kmax为zi=±z0时滞回曲线的最大斜率;Ffmax_z为滞回曲线在zi=±z0时两切线垂向距离的一半,如图7 所示[9]。库仑摩擦模型中的初始力Ffs_z和初始位移zs可由加载特性获取。
主气室静态试验低频加载时节流孔中气体流动速度缓慢,阻尼效果不明显[11],故结合式(6),已知空簧内压时气动刚度Kza为:
式中Ve(zi)为空簧有效体积。
则分数导数Kelvin-Voigt 模型中线弹性刚度Kze为:
主气室动态试验时胶囊黏弹性不可忽略,记动态特性试验时主气室动刚度为KD_Z,则分数导数模型中分数导数项对主气室总体刚度的贡献KK_V为:
式中Kzf为摩擦模型刚度与图7 中的“Kze+Kza”类似,KD_Z可由动态试验滞回曲线zi=±z0两点处斜率求得。已知KK_V时,式(3)分数导数项可等效表示为:
动特性加载曲线为正弦型,即
式中AD为加载位移幅值;ω为圆频率。式(38)可写为:
在求得不同频率对应的KK_V后,将其用形如下式的曲线拟合:
则根据下式可求得分数导数项的值:
应急止挡、胶囊材料、气动作用力所对应的刚度和阻尼按照式(31)和(32)计算。在已知滞回曲线时,式(32)中ΔE可表示为滞回环面积的函数,按文献[11]计算。
本文所建立的空气弹簧模型参数见表1。表1中,尺寸参数各符号与图3(a)对应,特征参数各符号与式(6)~(9)对应。需要说明的是,部分模型参数及后续实验数据源自文献[14-15],有效面积、有效体积及其变化率等部分特征参数通过几何关系计算得到,限于文章篇幅,此处仅给出计算结果。
表1 空气弹簧模型参数Tab.1 Parameters of air spring model
空气弹簧的垂向静、动态特性测试依据标准文献[16]执行,加载曲线如图8 所示。
图8 空气弹簧静、动态测试加载过程Fig.8 Loading process of static and dynamic test of air spring
静态测试时安装高度为标准工作高度260±2 mm,加载波形为三角波,加载速度为0.8 mm/s,采用10 mm 振幅时的结果评定空簧垂向刚度,在气囊内分别充入0.7,0.6,0.5,0.4,0.3 MPa 压缩空气,每个工况循环3 次,取第3 次测得的力-位移曲线,静刚度测试加载过程如图8(a)所示。动态测试安装高度为260±2 mm,加载幅值为70 mm 的正弦波,在气囊内分别充入0.6 和0.3 MPa 压缩空气,每个工况循环5 次,取第5 次测得的力-位移曲线,动刚度测试加载过程如图8(b)所示。
空气弹簧满载载荷Fm=13 kN,在进行应急止挡的力-位移特性测试时,最大加载载荷取1.25Fm[16],即16.25 kN。测试时,将0~16.25 kN 的载荷以10 mm/min 进行加载,两相邻载荷之间测试不停顿,进行3 次循环加载,从第3 次循环开始记录数据。
内压为0.3 MPa 时,空气弹簧主气室静态测试得到的力-位移滞回曲线如图9 所示。
图9 主气室静态测试滞回曲线Fig.9 Hysteresis curve of static test of main gas chamber
将图9 中的结果代入式(34)~(36),得:z2=1.22 mm;Kza=45.68 N/mm;Kze=3.91 N/mm;Kzf=9.95 N/mm。再将上述参数代入式(37)~(44),结合动态测试结果辨识分数导数项参数,得:az=0.1997,bz=2.1985。
应急止挡载荷-位移曲线如图10 所示。对其线性拟合,95%置信区间内其刚度可用斜率表示为:
图10 应急止挡载荷-位移曲线Fig.10 Emergency stop load-displacement curve
结合式(32),其阻尼为:
至此,垂向模型各参数辨识完毕。将各参数代入式(1)~(6),绘制出内压为0.3 MPa 时的各作用力如图11(a)所示;内压为0.3 和0.6 MPa 时垂向力Fz的计算值与测试值对比如图11(b)所示。图11(a),(b)表明,±70 mm 范围内,气动作用力Fza为空簧垂向力Fz的主要成分;大于70 mm 时,应急止挡力Fzs发挥主要作用。计算值与测试值能够较好吻合,验证了本文垂向模型及其参数辨识结果的正确性。但式(1)结构成分的复杂性导致其应用难度大,故拟通过拟合Fz-z曲线以建立图6(c)模型对应的函数表达式。
图11 空气弹簧测试结果及计算结果Fig.11 Test results and calculation results of air spring
图10 和11 表明,大于70 mm 时止挡力Fzs可线性表示,故以70 mm 为临界点对Fz进行拟合。zi∈[-70,70] mm 时,由式(1)可知,Fz由三项组成,考虑常数项,采用三次多项式拟合Fz-z曲线,0.3 和0.6 MPa 时拟合结果分别为:
其中,Fz单位为N,拟合度R2分别为0.9982和0.9987。
图11(c)为标准高度下不同内压对应的空簧载荷曲线。对测试值进行线性拟合,拟合度高达0.99987,这表明Fz-p呈线性关系。故可对式(47)进行线性插值,得到不同内压与±70 mm 行程范围内的垂向力为:
依据式(24),(25)和(27)计算得 到内压0.3 MPa 时各横向力如图12 所示。图12 中Fya和Fye分别为与胶囊变形和加载频率相关的力,Fys和Fyf几乎不受加载频率影响,当横向位移达到对应值后即会作用。
图12 内压0.3 MPa 时的横向力Fig.12 Transverse force at 0.3 MPa internal pressure
依据前述结果,根据式(31)和(32)计算出空气弹簧系统的刚度和阻尼特性如图13 所示。
图13 空气弹簧系统刚度、阻尼特性曲线Fig.13 Stiffness and damping characteristic curves of air spring system
图13(b)中,横向位移0~0.5 mm 为空簧上、下盖板最大横向相对位移,在此范围内主要为胶囊变形产生的气动刚度、摩擦刚度以及胶囊自身的黏弹性力对应的刚度;0.5~6 mm 为应急止挡变形区,刚度即为图10 中的拟合斜率b;大于6 mm 时系统横向刚度主要由刚度更大的中央导柱提供。
结合线路动态测试结果,基于UM 建立磁浮车辆动力学模型,对比本文模型与线性模型的差异。线性模型中的刚度与阻尼通过空气弹簧动刚度试验确定,即刚度由橡胶气囊和附加气室等效计算获得,阻尼由节流孔阻尼效应等效计算得到[8]。线性模型参数取值为:纵/横向刚度为2.6×106N/m,阻尼为1860 N·s/m;垂向刚度为5.3×104N/m,阻尼为3110 N·s/m。
轨道不平顺幅值为±2.5 mm,依据文献[22-23]反演。以空簧上部固接于车体的导轨安装座评估其隔振性能,10~80 km/h 运行速度下导轨安装座处各向加速度有效值(RMS)如图14(a)~(c)所示。由图14(a)~(c)可见,速度低于30 km/h 时本文模型与线性模型计算结果相当,大于30 km/h 时线性模型计算结果逐渐偏大,且速度大于60 km/h 时偏大程度进一步增大。由前文分析结果可知,空簧刚度、阻尼与其动位移相关,以垂向为例,取30 和60 km/h 速度下30 s 内导轨安装座处垂向动位移计算结果进行统计分析,数据采样总数6000 个,不同位移出现频次如图14(d)所示。对比图14(d)可见,速度为60 km/h 时两模型计算结果差异显著,主要是动位移大于平均数时线性模型对应数据点个数大于非线性模型。
图14 直线线路模型对比Fig.14 Comparison of straight line model
中低速磁浮车辆可通过最小曲线半径50 m,目前实际运营线路最小曲线半径100 m,全线最高运营速度100 km/h,为兼顾120 km/h 提速运营的需求,进一步考虑选取800 m 半径曲线,各曲线线路参数依据标准文献[23]设置。半径50,100,200,300,800 m 曲线对应最大通过速度15,40,60,70,120 km/h[23]。由于空簧纵/横向模型相同,且平曲线线路横向动态作用显著,故以下仅给出各曲线上导轨安装座处横向、垂向计算结果,如图15 所示。
图15 曲线线路模型对比Fig.15 Comparison of curve line model
R=100 m 时,速度为5,15,25 km/h 时的测试结果分别为横向0.05,0.18,0.56 m/s2,垂向0.02,0.07,0.12 m/s2;R=300 m 时,速度为20,40,60 km/h 时的测试结果分别为横向0.33,0.96,1.72 m/s2,垂向0.09,0.27,0.51 m/s2,与 图15(a)~(b)对比可见其与本文模型较为接近。由图15(a)~(b)可见,同一速度下(低速)曲线半径越小,线性模型计算结果偏差越大,同一曲线半径车辆通过速度越大,线性模型偏差越大。以车辆速度为60 km/h 通过200 m 半径曲线为例,如图15(c)所示(图15(c)中:区间①为直线;②为缓和曲线;③为圆曲线),可见导轨安装座处垂向加速度线性模型几乎全程均略大于本文模型,圆曲线处尤为显著;横向加速度线性模型在缓和曲线与圆曲线过渡处以及圆曲线处显著大于本文模型,直线处两者较为接近。这表明,在提速时车辆动力学性能预测以及100 m 以下小半径曲线工程应用中线性模型的计算结果存在较大误差。
进一步对比两种模型的频率特性。车辆以15 km/h 通过R105 m+45‰曲线线路时,经空簧减振后导轨安装座处各向振动加速度功率谱密度如图16 所示。
图16 导轨安装座处各向振动加速度功率谱密度(PSD)Fig.16 Power spectral density(PSD)of vibration acceleration in each direction at guide rail mounting seat
图16(a)为0~1000 Hz 内导轨安装座处的垂向PSD 图,可见低频内两模型PSD 幅值和变化趋势均存在显著差异,但随着频率增加两者变化趋势逐渐接近,线性模型PSD 值略大于非线性模型。空簧系统固有频率一般在0.7~3 Hz[21],该型空簧固有频率f0=2.62 Hz,当激励频率大于时发挥隔振效果,故选取0~6 Hz 对比分析两模型频率特性的差异。由图16(b)~(d)可见,纵/垂向频率小于1 Hz、横向频率小于5 Hz 时线性模型PSD 值显著小于非线性模型。结合前文空簧系统刚度、阻尼特性曲线,可认为两模型出现差异的原因是:线性模型恒定的刚度与阻尼导致其对低频减振特性表征不足而过高地评价了中高频减振特性。
基于振动力学与弹性力学基本原理,构建了考虑胶囊、中央导柱、应急止挡等关键部件力学特性的中低速磁浮车辆空气弹簧垂向、横/纵向动力学模型,推导了垂向-横/纵向力学关系,探究了系统参数的辨识方法,基于力学等效原理给出了简化模型,并以某型空气弹簧为例,通过空簧试验与车辆线路动态测试结果对比验证了本文模型的准确性与优越性,主要结论如下:
(1)建立了考虑胶囊黏弹性特性、中央导柱摩擦/支撑特性以及整体气动特性的磁浮车辆空气弹簧垂向、横/纵向非线性模型,构建了垂向与横/纵向的力学关系,探究了系统参数的辨识方法,试验验证了模型的准确性。
(2)±70 mm 行程范围内,磁浮车辆空气弹簧垂向力与内压、高度之间的关系可用三次多项式准确表述;垂向位移大于70 mm 时,载荷与位移呈线性关系,结构的特殊性导致磁浮车辆空气弹簧横向刚度远大于垂向,且近似呈分段线性。
(3)磁浮车辆提速时的动力学性能预测以及半径100 m 以下曲线线路动态通过性能计算时线性模型误差较大。具体表现为:直线线路上车辆速度大于30 km/h 时线性模型计算结果大于本文模型,且速度越高偏差越大;半径小于100 m 的曲线线路上线性模型计算结果显著大于本文模型;同一半径曲线线路车辆通过速度越大线性模型偏差越大,车辆速度大于60 km/h 时尤为显著。刚度与阻尼恒定的线性模型对低频减振特性的表征不足和对高频减振特性的过高估计是造成计算误差的主要原因。