车线动力学(耦合振动3)

1、 第五节半车第五节半车轨道空间耦合轨道空间耦合 振动模型的有限单元解法振动模型的有限单元解法 n车辆轨道空间耦合振动求解，涉及到 系统中众多的： 耦合关系 非线性因素 时变特性 n求解过程十分复杂，因而最有效的方法 为有限单元法。 一、模型特点及变量选择一、模型特点及变量选择 半车轨道的“垂向横向非线性空间耦合振动时变模 型” Mc,Jc n耦合特性： 轮轨垂向振动耦合：轮轨接触力或弹簧 轮轨横向振动耦合：蠕滑力和轮缘力 n非线性特性(坐标原点位于系统不受任何力时 的质心位置)： 轮轨间垂向力或非线性赫兹接触弹簧 轮轨间非线性蠕滑力 n时变特性： 半车模型在轨道移动，车辆每前进一步， 振动方程

2、均需重新组建。 车辆自由度 n车体：沉浮(车体在曲线上受到的未被平衡离 心力则以等效离心力和力矩施加于转向架二系 悬挂点上) n转向架构架：横移、沉浮、摇头、侧滚及点头 n轮对：横移、沉浮、摇头、侧滚 n车辆自由度共14个。 轨道自由度 n钢轨：8N(每一枕跨为一个单元，每一钢轨节 点有垂向位移、绕y轴的转角、横向位移、绕z 轴的转角4个未知变量，两股钢轨) n轨枕：3N(沉浮、横移、侧滚) n轨道自由度数：11N n轮轨间设置4个轮缘力未知变量。 (1)当轮缘贴靠钢轨时，轮缘力未知量列 入方程，同时引入轮轨横向位移协调条件； (2)而当轮轨不贴靠时，总刚中与轮缘力 对应的主元赋大值。 n振动

3、模型中总的未知变量数为11N+18个。 变量位置 部件位置 曲线外轨(4N)14N 曲线内轨(4N)4N+18N 轨枕(3N)8N+111N 第1轮对(4)11N+111N+4 第2轮对(4)11N+511N+8 构架(5)11N+911N+13 车体(1)11N+14 轮缘力(4)11N+1511N+18 二、系统动能、势能、虚功二、系统动能、势能、虚功 及振动方程的组建及振动方程的组建 、振动系统的总动能、振动系统的总动能( (质量矩阵质量矩阵) ) e i e N Ni i eT i e i e N i i eT iumuumuT . 12 1 . 1 1 . 2/12/1 左右钢轨 2

4、/1 2 1 2 2 si s N i si s si s Jmm y z 轨枕 JJmm1/2 2 wi w 2 wi w 2 wy w 2 wi 2 1i wzy 轮对 z y zc c b b b b b b b b b b mJJJmm 2 2222 2 2/1 构架车体 . . 振动系统的总应变能振动系统的总应变能( (刚度矩阵与荷载列阵刚度矩阵与荷载列阵) ) （未包含轮轨接触弹簧的应变能） e i N Ni e i eT i e i e i N i eT i ukuukuU 12 1 1 1 2/12/1 外轨弯曲内轨弯曲 )()(2/1 2 )( 2 1 1iNrsirisi

5、N i r yyyyk 外轨扣件 横向刚度 内轨扣件 横向刚度 N i iNrsisirisisir zbzzbzk 1 2 )( 2 2 )()2/1 外轨扣件 垂向刚度 内轨扣件 垂向刚度 )()(22/1 2 2 2 2 1 2 1sisissisis N i si s bzkbzkk y 道床横向 刚度 道床外侧 垂向弹簧 刚度 道床内侧 垂向弹簧 刚度 4 1 2 2 1 2 1 41 2 1 1 2 ) 1(2/1 i i w bb i by y R Hlyk l 一系悬挂横向刚度 4 1 2 1 1 112 13 1 2 1 1 1 )1()1()1(2/1 i b i i i

6、w i x b R lb bk 一系悬挂纵向刚度 曲线半径的 影响，进入 荷载列阵 2 2 11 4 1 2 111 2 1 1 ) 1() 1() 1(2/ 1 i w i i i w b i b i bz bzblzk 一系悬挂垂向刚度 2 2 2 2 2 2 2/12/12/1 bbb kkk 二系悬挂摇 头角刚度 二系悬挂侧 滚角刚度 二系悬挂点 头角刚度 yb ybcz kzzk 2 2 2 2 2/1)(2/1 二系悬挂 垂向刚度 二系悬挂 横向刚度 、振动系统中的阻尼力虚功、振动系统中的阻尼力虚功( (阻尼矩阵阻尼矩阵) ) )( )( )()( 1 1 1 11 iNrsiiN

7、rsi N i r risirisi N i r yyyyc YYyycW 外轨扣件 横向阻尼 内轨扣件 横向阻尼 )( )( )( )( )( 1 2 1 2 iNrsisi iNr si si N i r risisi N i ri si si r zbzzbzc zbzzbzc 外轨扣件 垂向阻尼 内轨扣件 垂向阻尼 )( )(2 1 2 2 1 1 sisisi si N i s sisisi si ssi N i sis bzbzc bzbzcYyc 道床横向 阻尼 道床外侧 垂向阻尼 道床内侧 垂向阻尼 ) 1() 1() 1( ) 1() 1() 1( 2 11 2 111 2

8、1 2 1 1 4 1 2 1 11 2 1 1 i w i i w b i b i b i w i i i w b i b i b z bzblz bzblzc 一系悬挂 垂向阻尼 bbybbbbbb bc bc z yycccc zzzzc 2 . 2 . 2 . 2 2 )( 二系悬挂 垂向阻尼 二系悬挂 摇头阻尼 二系悬挂 侧滚阻尼 二系悬挂 点头阻尼 二系悬挂 横向阻尼 4 4、蠕滑力虚功、蠕滑力虚功( (荷载列阵荷载列阵) ) 4 1 2 11 1 2 ) 1( i i w i i bTW 纵向蠕滑力 2 10 4 1 2 12 i w ri i i w i rYYT横向蠕滑力 2

9、 13 4 1 i w i i M自旋蠕滑力矩 5 5、轮缘力及轮缘摩擦力虚功、轮缘力及轮缘摩擦力虚功( (刚度矩阵刚度矩阵) ) 4 1 4 1 2 11 1 3 ) 1() 1( ii i w i ie iii i YFuNFW 轮缘横向力 2 1 4 1 1 ) 1( i w i i i bF 4 1 2 10 1 ) 1( i i w i i rF 轮缘摩擦力 6 6、曲线上未被平衡的离心力所作的虚功、曲线上未被平衡的离心力所作的虚功 ( (荷载列阵荷载列阵) ) 2 1 314 )(5 . 0 i bbcbbwiwi HHyFyFyFW 轮对构架 车体 7 7、轮轨接触力虚功、轮轨接

10、触力虚功(刚度矩阵及荷载列阵) 2 2 1 4 1 2 1 ) 1(2/1 ririi w i i i w cic zbzkU 处理成弹簧的应变能 ) 1( 4 1 2 1 2 15 i i w i i w iriirii bzPPzPW 处理成力的虚功 垂向力 ri i yii w i yi yPyP 4 1 2 1 4 1 4 1 2 10 4 1 2 1 2 1 i i w yi i i w i w yi rPbP 垂向力的 水平分力 8 8、重力所作的虚功、重力所作的虚功(荷载列阵) 2 1 6 5 .0 i wiwbbcc zgmzgmzgmW 三、振动方程的求解方法及过程三、振动方

11、程的求解方法及过程 n振动系统的非线性振动方程组： ),(),()()( uutPutFtuCtuM 不同时刻车轮处于轨道的不同位置，钢轨不平 顺值也不相同，荷载列阵发生变化。 荷载列阵中包含蠕滑力，隐含振动位移、速度。 轮轨接触力是波动的，轮轨接触弹簧瞬时刚度 依赖于接触力。 n循环顺序 时间及车轮位置 轮缘与钢轨的贴靠判断 接触力迭代 蠕滑力饱和判断及蠕滑力迭代 1 1、选择积分方法、选择积分方法 n根据方程组的强非线性和时变特性，选 取无条件稳定且无超调现象、对非线性 方程有很强适应性的Park积分方法，对 微分方程进行求解。 nPARK方法的表达式： 3 1 11 0 1 3 1 1

12、1 0 1 i in i nn i inin n uauau uauau t a t a t a t a 6 1 , 1 6 15 , 6 10 32 10 n代入振动方程可得到： ),()( )()( 1 . 11010 . 0001 2 0 n nnnn n nn uuPFuBuaC uBuBauM a 3 1 1 0 3 1 10 )( ,)( i in i n i inin uauBuauB 、处理非线性因素、处理非线性因素 n振动系统中主要包含非线性接触刚度和 非线性蠕滑力两种非线性因素。根据这 两种非线性因素不同的特点，采用不同 的迭代方法进行处理。 （1 1） 用直接迭代法处理非

13、线性蠕滑力用直接迭代法处理非线性蠕滑力 n因蠕滑力与 蠕滑的关系 曲线是凸曲 线，蠕滑力 的直接迭代 过程是稳定 和收敛的。 0 n对于每一次蠕滑力迭代，振动方程中的 荷载列阵变为常量。 )( )()( 1010 . 0001 2 0 PFuBuaC uBuBauM nnn n nna （2 2） 用牛顿迭代法处理非线性接触力用牛顿迭代法处理非线性接触力 n轮轨非线性接触 力与位移的关系 曲线是凹曲线， 直接迭代方法是 发散的。因迭代 中非线性刚度取 为切线刚度，改 善迭代收敛性的 最简化途径是采 用牛顿迭代法。 ) ( ) ( 1 )1 ( 1 ) ( 1 u u K i n i n i n

14、T 0)( )()( )()( 0 . 0001 1101 2 0 1 n n nn nnnn uBC uBMuBMaP uFuCauMu a ) ( 1 )1 ( 1 ) ( 1 uuu i n i n i n )( 1 )( 1 )( 1 n i n i nT u u K )( 1 )( 1 0 2 0 n i n u uF CaMa 、启步方法、启步方法 nPark方法是线性多步法，不能自动启步， 需借助其它方法启步。自动启步的方法 很多，可选取最常用的Wilson法。 n启步将轮轨接触简化为线性弹簧，可考 虑非线性蠕滑力。使启步运算比较简便， 避免对于非线性方程不收敛。 )( 6 5

15、. 0 )( 2 )( 1 2 2 1 1 1 nnn nn nnnn nnn uu t ututuu uu t utuu uuuu n振动方程组变为： )()()(PtuKtuCtuM BuA 3 1 tC(0.5- 6 5.0 2 2 nn nn n uKtu utuCuKPB K t CtMA 4 4、程序编制、程序编制 n计算程序包含一个主程序和多个子程序 （1）钢轨有限和无限单元刚度矩阵； （2）钢轨有限和无限单元质量矩阵； （3）钢轨梁无限单元阻尼矩阵； （4）组集总质量矩阵； （5）计算各车轮位置、车轮所在钢轨 单元形函数值、车轮下轨面不平顺值； （6）处理各部弹簧及相关的常数项

16、； （7）处理各部阻尼及相关的常数项； （8）处理曲线上求被平衡的离心力、重力、 轮轨间垂向力的水平分力、轨面不平顺； （9）计算轮轨间的蠕滑力。 （10）判断轮缘与钢轨的贴靠状态，处理 轮缘力及轮缘摩擦力； （11）进行蠕滑力饱和判断及修正； （12）用全主元高斯方法求解代数方程； （13）时间为0时计算各单元初始位置； （14）第一、二步Wilson方法启步； （15）从第三步开始用PARK方法处理振动 方程。 n程序外层 结构框图 t=t,2t Wilson-方法启步 t=0 求解静平衡位置 开始 t=3t.nt Park方法解振动方程 结束 组建总质量和 总阻尼阵 计算单元矩阵 KP=

17、1 B:KP=1 轮轨接触力取切线刚度 KP=2 轮轨接触力为常值 D:用Park方法处理 振动方程 E:用高斯法解方程 F:KP=1 判断位移增量 平方和小于 给定值 否 是 G:判 断轮缘贴 靠状态 否 是 H:KP=1 轮轨接触力计算 C:组建轮轨系统 振动方程 A:计算各车轮位置及轮轨接 触点的钢轨单元形函数值 I:KP=1 判断接触力增 量平方和小于 给定值 KP=2 K:判断蠕滑饱 和及修正蠕滑力 否 是 否 是 L:判断 蠕滑力平方 和小于给 定值 J:计算蠕滑力 输出结果 M:计算磨耗功 nPARK方 法程序 框图 程序功能 n程序中线路曲线半径、轨道垂向和横向 不平顺可任意选

18、择，且内外轨不平顺可 以不一样。 n运算结果输出为轮轨间蠕滑功、蠕滑力、 轮轨垂向力、轮缘力，轨道振动的力、 位移、速度和振动加速度。 5 5、计算方法的参数选择、计算方法的参数选择 n正式进行计算之前，除对车辆、轨道系 统中各有关结构参数要进行研究和选择 外， n须对对程序计算结果进行校验， n并进行数值试验，选择与计算方法有关 的合适的计算参数，如积分步长、收敛 精度、模型长度等。 （1 1）积分步长）积分步长 nPark方法是无条件稳定的，任意积分步 长都能计算出结果。 n但步长增加时算法阻尼增大，步长太小 耗费机时且增加舍入误差，因此必须根 据求解的频率范围选择合适的步长。 n根据Pa

19、rk方法算法阻尼与积分步长的关 系，对应振幅误差率小于2.5的最高振 型频率为0.1/t（Hz），对应振幅误差 率小于15的最高振型频率为0.15/t （Hz）。 n如积分步长取为0.001秒，则对应100Hz 振型的振幅误差率为2.5，对应150Hz 的振型的振幅误率为15。 n如要确保400Hz以内的振幅误差率小于 2.5，则积分步长应小于0.00025秒。 n实际计算中，为节省机时，根据不平顺 的长度和列车速度进行时间步长数值试 验。步长从0.001秒开始反复折半，计算 出23组结果，当前后两组结果的误差 在510以内时，认为步长是合适的。 （2 2）迭代误差限）迭代误差限 n迭代误差限