基于改進的Newmark-β法重載列車車鉤縱向力仿真研究

武承龍，董昱

基于改進的Newmark-法重載列車車鉤縱向力仿真研究

武承龍，董昱

(蘭州交通大學 自動化與電氣工程學院，甘肅 蘭州 730070)

重載列車運行過程中過大的車鉤縱向力一直是制約重載列車發展的瓶頸，空氣制動不同步是產生列車縱向沖動的根源，導致車體擠壓車鉤形成車鉤力。傳統的經過制動特性試驗采集車鉤力的方法耗時耗力，為了經濟地獲取重載列車在不同線路上運行時車鉤力的大小，將Newmark-法應用于重載列車車鉤縱向力的仿真分析中。由于列車縱向動力學方程是非常復雜的非線性方程，傳統方法為了保證計算精度而采用大量迭代運算，耗時長效率低。基于增量思想改進Newmark-法，通過引入預測解直接對非線性方程進行處理，然后對預測解進行校正，最終得到收斂的近似解。算例結果表明，改進算法在保證了計算精度的同時計算效率更高，更適用于長大編組重載列車車鉤縱向力的仿真計算和分析。

重載列車；Newmark-法；車鉤緩沖器；縱向力；預測解

重載運輸是提高運輸效率、降低運營成本、提高鐵路綜合競爭力和爭取市場份額最為有效的方式之一。我國鐵路從20世紀80年代開始發展重載運輸，至今已建成包含神華集團的神朔鐵路、朔黃鐵路以及大秦線在內的一系列重載鐵路運煤專線。“貨運重載”成為繼“客運高速”后，中國鐵路建設的重點方向和爭取市場份額最為有效的方式之一。然而由于重載列車空氣制動系統自身的局限，編組越長的列車制動波傳遞速度也就越慢，就會造成列車前后部制動的不同步，后車擠壓碰撞前車產生過大的車鉤縱向力，導致車鉤斷鉤的重大安全事故。重載鐵路線路多修建于惡劣山區路段，線路條件的復雜也造成了重載列車受力情況復雜，難以建立數學模型對車鉤縱向力精確計算。為了進一步分析研究重載列車車鉤縱向力對列車運行過程的影響，Duncan等[1]將緩沖器加載曲線簡化為兩段線性模型，由于模型較為簡單，難以從數學角度精確描述。孫樹磊等[2]構建了適用于貨車摩擦緩沖器的多段線性動力學修正模型，利用車輛沖擊試驗驗證了模型的正確性。魏偉等[3]開發了基于空氣制動系統仿真的列車縱向動力學仿真平臺，計算了2萬t列車車鉤力的分布特性；由于列車的運行過程均呈現出強非線性特性，Newmark-法也開始被廣泛應用于列車運行問題的求解。Fancher等[4]采用Newmark顯式積分方法求解非線性的列車動力學方程；ZHANG等[5]提出了一種精細積分法，用于非線性輪軌接觸的車輛?軌道耦合動力學分析中，具有很高的計算精度和準確性；張斌等[6]通過對車輛?軌道垂向耦合系統動力響應計算的迭代過程進行了改進，在一定程度上提高運算效率；常崇義等[7]提出了基于Newmark-的高精度平衡迭代算法，通過迭代運算提高了算法的精度，但是計算過程復雜，時效性不能得到保證；本文通過對重載列車運行過程中的縱向受力情況進行分析，結合車鉤緩沖器的工作原理及其多段線性特性，得出緩沖器加載曲線函數和卸載曲線函數。為了避免對于重載列車車鉤縱向力計算過程中非線性項的反復線性化過程，本文提出了一種基于增量思想的預估校正的計算方法，增量的變化方式采用二分法，目的是快速逼近收斂解。通過計算數據和實驗數據比對，驗證了改進算法的有效性，保證了計算精度的同時，提高了計算效率。

1 重載貨車車鉤緩沖器動力學模型

1.1 鋼彈簧與干摩擦組合式緩沖器結構

鋼彈簧與干摩擦組合式緩沖器是我國鐵路上應用最為廣泛的緩沖器之一，這種緩沖器構造簡單，成本較低并且緩沖效果好。國內代表性的該類緩沖器有MT-2型車鉤緩沖器和MT-3型車鉤緩沖器，其構造如圖1所示。加載過程中，車鉤連接處推動楔塊左右移動，楔塊上部分緊壓斜板通過摩擦消耗一部分沖擊產生的能量，還有一部分能量轉化為彈簧的彈性勢能。在卸載過程中，由于彈簧處于壓縮的狀態下，楔塊被推回到原位置，楔塊與斜板之間仍可通過摩擦消耗列車縱向沖擊產生的能量，這樣來回的加載和卸載過程就可將部分縱向沖擊力消耗。

圖1 彈簧與干摩擦組合式緩沖器

1.2 車鉤緩沖器多段線性模型

由以上分析可知，在緩沖器加載和卸載的過程中，楔塊同斜板的接觸面積、摩擦角度及正壓力在時刻改變，因此精確構建表示緩沖器特性曲線的數學模型有非常大的難度。緩沖器制造廠商均會提供車鉤緩沖器落錘試驗的動態特性曲線。可以獲得車鉤緩沖器在不同行程處返回平衡位置的動態特性曲線及對應加載狀態和卸載狀態的動態包絡線。如圖2所示為MT-2型車鉤摩擦緩沖器的落錘試驗特性曲線[8]。

圖2 MT-2型緩沖器落錘試驗及沖擊試驗曲線

由于在列車運行的過程中，車鉤緩沖器同時存在有彈簧的彈性力及阻尼力，彈簧力與相對位移有關，而阻尼力與相對速度大小成正比，但與相對速度相反。因此可得到緩沖器的阻抗力如式(1)所示[9]。

其中：

圖3 MT-2型車鉤緩沖器多段線性模型

2 改進的Newmark-β法對重載列車縱向力學方程的求解

2.1 Newmark-β法原理與步驟

對于縱向力學方程[10]：

Newmark-法根據時間步長內加速度的變化規律來對最終輸出響應進行計算，通過對時間步長內加速度進行假設，以時刻的運動量為初始值，通過積分計算得到一個時間步長+Δ的運動公式。

加速度同時由2個參數調整和控制：

+Δ時刻的速度和位移：

聯立以上2個方程組，可得+Δ時刻的速度和加速度的計算公式：

由于運動仍需滿足+Δ時刻的縱向動力學平衡方程：

由式(8)，式(9)和式(10)可得：

將式(11)寫成如下形式：

求解步驟如下：

2) 確定系統的有效剛度矩陣；

3) 對+Δ時刻的有效載荷進行求解；

4) 用有效剛度矩陣和有效載荷矩陣對+Δ時刻的位移求解：

5) 對+Δ時刻的速度和加速度進行求解：

6) 重復以上步驟，可計算下一個時間步長Δ的運動情況。

2.2 系數矩陣的確定

1) 剛度矩陣[]

剛度矩陣的概念源于固體力學，是利用有限元的方法計算單元受力的一個重要系數矩陣[11]。剛度是表示物體抵抗形變的能力，其物理含義是當某個組成整體的單元位移發生改變時，產生的作用于節點處力的大小。如圖4所示，是一個3節車體之間通過車鉤緩沖器相連接組成的體系，緩沖器的剛度為，車體位移u，受力f。

圖4 重載列車剛度矩陣的構建

根據胡克定律可得：

寫為一般的形式：

其中：[]為該體系的整體剛度矩陣，其每一個元素k表示第個單元在坐標軸方向上發生單位位移后，第個單元保持零狀態時所需要施加的力的大小。k的取值根據式(2)，式(3)和式(4)確定，首先計算相鄰兩節列車的相對位移差，確定緩沖器處于哪種工作狀態，例如：當緩沖器處于加載狀態的時候，根據前后位移差選擇所處的線性段，選取式(2)中對應的k(=1,2,…,)值。當系統由個單元組成時，剛度矩陣[]為×的方陣。

2) 阻尼矩陣[]

阻尼矩陣[]是由質量矩陣[]和剛度矩陣[]成的。目前在工程應用中，阻尼矩陣常用Rayleigh阻尼陣的形式來構造，即：

3) 質量矩陣[]

質量矩陣即各單元質量構成的矩陣。對于重載列車車鉤縱向力的計算，質量矩陣由每一節貨物列車滿載時的總重和機車自重組成。

2.3 改進的Newmark-β法

本文采用直接對非線性方程進行處理的方法，跳過將非線性問題轉化為近似的線性問題這一個步驟，省略了復雜的反復迭代運算過程，大大降低了計算量和累計誤差的產生。

設時刻系統的狀態量為：

經過時間Δ后系統的狀態量為：

針對于+Δ時刻有效載荷無法精確計算，用以下步驟進行改進：

2)在時刻系統{}的基礎上加上增量Δ，得+Δ時刻系統的猜測解：

5) 定義截斷誤差的表達式：

將代入上式計算方程的殘差Re即不平衡力，若不平衡力大于允許的限值且為正，則返回第一步，利用二分法重新設定Δk，另重新對t+Δt時刻的有效載荷進行預測，再次計算不平衡力，若不平衡力時，則t+Δt時刻的解，然后進入下一個時間步的計算，本文算法流程圖如圖5所示。

2.4 算法的收斂性和穩定性

Newmark-法的控制參數和的取值，與算法自身的精度和穩定性有著重要的聯系，當=1/2時，方程才具有二階精度，因此一般均取值：

系統穩定的條件：

式中：Δ為時間步長；T為振動系統的最小固有 周期。

3 重載列車緊急制動過程運動狀態及車鉤力的計算分析

本文重載列車編組方式采用“1+1”編組，即“HXD1型電力機車+105輛C80型貨物列車+HXD1型電力機車+105輛C80型貨物列車+可控列尾裝置”。機車選用HXD1型電力機車，制動系統采用CCBⅡ型空氣制動系統，整備質量2×92 t，裝有13號車鉤及MT-3型緩沖器。貨物列車選用C80型煤礦專用敞車，載重80 t，自重≤20 t(按20 t計算)，空氣制動裝置采用120型控制閥，采用高摩合成閘瓦。車輛1位端安裝E級鋼17號車鉤，采用MT-2型緩沖器。

試驗條件及運行工況為“1+1”編組滿載重載列車以80 km/h的速度運行至大秦線k142+216處(坡道坡度：8‰，下坡道)進行緊急制動。本文以下仿真計算和分析均是基于此工況及相關參數展開。

3.1 重載列車緊急制動過程中加速度和速度的數值計算和分析

制動特性是導致列車縱向沖動產生的根本原因。列車不同部位開始減壓時間不同，制動波的傳遞速度也就不同。圖6所示是“1+1”編組重載列車制動缸的動作時序。列車首中尾部各對應一個最小值，機車位置及可控列尾裝置安裝處附近車輛制動缸減壓速度要比中部機車速度快，制動波傳遞也是從機車所在部位向中部貨物列車依次傳遞。圖7和圖8記錄了從22 s處開始施加緊急制動命令，制動閘瓦壓緊車輪開始，到整列車完全停住的過程中，第1，52，107，172和212輛車緊急制動的速度曲線和加速度變化曲線。

由圖7和圖8及運算數據可得列車緊急制動過程中，制動時間38.6 s，制動距離498.27 m，單節車體減速度范圍在0.473～0.886 m/s2之間。根據《牽引計算規程》中提供的計算公式及方法，對算例中的重載列車制動過程進行計算，牽引計算的結果為：以80 km/h的速度在坡度為8‰的下坡道上進行緊急制動，制動距離469 m，整列車制動減速度范圍0.531～0.584 m/s2之間。對比牽引計算和縱向力學方程解算的結果，制動距離存在一定的誤差，考慮到牽引計算是將整列車看作一個質點而進行的計算，忽略了列車之間復雜的相對運動和碰撞，只考慮了制動力和附加阻力的影響；而動力學方程解算的過程中，將每一節車運動方程單獨進行求解，由于車鉤縱向力的存在，作用于每一節車體的力變化明顯，單節車體減速度范圍波動較大。

圖6 制動缸動作時序

圖7 重載列車速度響應時程

實施緊急制動前，列車處于惰行工況，由于進行數值計算時初值的選取采用列車惰行的實測速度數據，實測數據表明每節車輛的瞬時速度都是不同的，相鄰車輛瞬時速度差在0.13～1.72 m/s的范圍內變化，經過單位換算呈現出如圖7所示的緊急制動前不同車輛間存在速度差的情況。前后車之間速度差的存在導致在一個時間步長Δ的范圍內位移也存在不同，車輛擠壓車鉤產生車鉤力作用于相鄰車輛，造成車輛受力的變化，出現圖8所示緊急制動前不同車輛之間存在加速度差。

圖8 重載列車加速度響應時程

緊急制動過程的前半段，第1位，107位和212位車加速度變化趨勢基本一致，加速度的大小在緊急制動起始時刻呈現一個先增大后減小的趨勢，是由于閘瓦摩擦因數是和速度相關的一個參數，速度越大閘瓦摩擦因數越小，制動力也就越小。隨著制動進行，空氣制動力越來越大，制動減速度也相應增大。先制動的機車速度已經降低，后制動的貨物列車速度大于機車，發生后車碰撞擠壓前車的現象，抵消一部分制動力，整車受力隨即減小，加速度減小。中部貨物列車，例如第52位和第172位列車，加速度大小變化均呈現一個越來越大的趨勢。當中部貨物列車碰撞前面已經開始制動的車輛后，前車對后車也存在一個反作用力，這個力方向和制動力的方向一致，阻礙列車的運行，整列車受力增大，加速度增大。

3.2 重載列車緊急制動過程中車鉤縱向力的計算驗證

為驗證改進的Newmark-法應用于重載列車車鉤縱向力分析的有效性和可行性，本文采用中國鐵道科學研究院機車車輛研究所發布的大秦線HXD1型機車牽引2萬t組合列車試驗報告中相關試驗數據進行比對。如圖9所示記錄了第107(從控機車)，124，157和172輛車車鉤縱向力的仿真計算結果。實線為大秦線“1+1”編組重載列車在k142+ 216(坡道坡度：8‰，下坡道)處進行制動試驗的實測數據。虛線為基于改進的Newmark-法對“1+1”編組重載列車進行仿真計算的結果。

圖9 車鉤縱向力數值積分計算結果

根據比對可知基于改進的Newmark-法計算結果和文獻試驗結果接近，2條曲線變化趨勢大致相同。因此，改進算法用來計算重載組合列車的車鉤縱向力具有足夠的精度，在一定程度上可以模擬重載列車制動過程，模擬分析在各種線路運行條件下列車所受車鉤縱向力的大小。

3.3 改進的Newmark-β法和Newmark-β高精度平衡迭代算法對比分析

圖10所示為2種方法仿真計算得到的重載列車緊急制動過程中最大車鉤縱向力的分布情況。經過對比分析和計算，2條曲線平均偏差率為5.916%，曲線基本吻合。改進的Newmark-法雖然簡化了計算步驟，但是計算精度上并沒有太大的損失。車鉤縱向力偏差最大為第31位車處，偏差為123.7 kN。出現這種誤差的原因是增量Δ的取值過大，不平衡力剛好滿足設定的收斂范圍，所以將此時的Δ看作是這一個時間步長內的收斂解。但是最大誤差123.7 kN相比于第31位貨物列車所受最大車鉤縱向力1 105.5 kN以及列車緊急制動所受的最大制動力7 152.94 kN來說，占比分別為11.189%和1.729%，誤差造成的影響是在可接受范圍內的。

圖10 最大車鉤縱向力隨車長分布曲線

圖11 計算時間對比

圖11所示為2種算法應用與車鉤縱向力計算過程中所耗費的計算時間對比，改進的算法可通過調整單位步長內的增量，快速逼近收斂解，從而減少反復迭代計算的次數。在“1+1”編組重載列車以80 km/h的速度進行緊急制動的仿真計算中，當程序循環迭代計算7 000次時，改進算法相比于傳統算法計算效率提升50.097%，而且隨著迭代次數的上升，計算效率越來越高。當處理長大列車編組、仿真模擬線路參數長距離運行以及需要得到更為精確的模擬參數作為參考的情況下，改進的算法更適用于這類場景。

4 結論

1) 基于增量思想改進了Newmark-算法，并對“1+1”編組重載列車在8‰下坡道緊急制動進行了仿真計算，結果證明改進算法應用于重載列車縱向動力學方程求解的可行性，在保證計算精度的同時提高了計算效率。

2) 將仿真計算得到的制動過程中車鉤縱向力和文獻試驗結果進行比對，結果表明改進算法的應用在一定程度上可以模擬分析重載列車在各種線路運行條件下列車所受車鉤縱向力的大小。

3) 對縱向動力學方程的求解，用預測解代替了線性化后再反復積分消除誤差的繁瑣步驟，省略了線性化的過程，為類似的非線性問題的處理提供了新的思路和方法。

[1] Duncan I B, Webb P A. The longitudinal behaviour of heavy haul trains using remote locomotives[C]// Proceedings of 4th International Heavy Haul Conference, Brisbane, 1989: 587?590.

[2] 孫樹磊, 丁軍君, 周張義, 等. 重載列車縱向沖動動力學分析及試驗研究[J]. 機械工程學報, 2017, 53(8): 138?146. SUN Shulei, DING Junjun, ZHOU Zhangyi, et al. Analysis and test of heavy haul train longitudinal impulse dynamics[J]. Journal of Mechanical Engineering, 2017, 53(8): 138?146.

[3] 魏偉, 趙連剛. 兩萬噸列車縱向動力學性能預測[J]. 大連交通大學學報, 2009, 30(2): 39?43. WEI Wei, ZHAO Liangang. Prediction of longitudinal dynamic coupler force of 20 000 ton connected train[J]. Journal of Dalian Jiaotong University, 2009, 30(2): 39? 43.

[4] Fancher P S, Ervin R D, Macadam C C, et al. Measurement and representation of the mechanical properties of truck leaf springs[C]// SAE Technical Paper Series. 400 Commonwealth Drive, Warrendale, PA, United States: SAE International, 1980.

[5] ZHANG J, GAO Q, TAN S J, et al. A precise integration method for solving coupled vehicle-track dynamics with nonlinear wheel-rail contact[J]. Journal of Sound and Vibration, 2012, 331(21): 4763?4773.

[6] 張斌, 羅雁云, 雷曉燕. 改進迭代過程的車軌耦合振動數值解法及應用[J]. 工程力學, 2016, 33(3): 128?134. ZHANG Bin, LUO Yanyun, LEI Xiaoyan. Numerical approach for vehicle-track coupling vibration based on improved iterative process and its application[J]. Engineering Mechanics, 2016, 33(3): 128?134.

[7] 常崇義, 王成國, 馬大煒, 等. 2萬t組合列車縱向力計算研究[J]. 鐵道學報, 2006, 28(2): 89?94. CHANG Chongyi, WANG Chengguo, MA Dawei, et al. Study on numerical analysis of longitudinal forces of the 20 000 t heavy haul[J]. Journal of the China Railway Society, 2006, 28(2): 89?94.

[8] 趙旭寶, 魏偉, 張軍, 等. 緩沖器分段阻抗特性對重載列車縱向沖動的影響[J]. 鐵道學報, 2017, 39(10): 33? 42. ZHAO Xubao, WEI Wei, ZHANG Jun, et al. Influence of segment impedance characteristics of draft gear on longitudinal impulse of heavy haul train[J]. Journal of the China Railway Society, 2017, 39(10): 33?42.

[9] 孫樹磊. 重載列車縱向沖動動力學研究[D]. 成都: 西南交通大學, 2014. SUN Shulei. Research on heavy haul train longitudinal impluse dynamics[D]. Chengdu: Southwest Jiaotong University, 2014.

[10] 翟婉明. 非線性結構動力分析的Newmark預測?校正積分模式[J]. 計算結構力學及其應用, 1990, 7(2): 51? 58. ZHAI Wanming. The predictor-corrector scheme based on the newmarks integration algorithm for nonlinear structural dynamic analyses[J]. Computational Structural Mechanics and Applications, 1990, 7(2): 51?58.

[11] 劉廣, 劉濟科, 陳衍茂. 基于Wilson-和Newmark-法的非線性動力學方程改進算法[J]. 計算力學學報, 2017, 34(4): 433?439. LIU Guang, LIU Jike, CHEN Yanmao. An improved algorithm for nonlinear dynamic systems based on Wilson-and Newmark-method[J]. Chinese Journal of Computational Mechanics, 2017, 34(4): 433?439.

[12] Piechowiak T. Pneumatic train brake simulation method[J]. Vehicle System Dynamics, 2009, 47(12): 1473?1492.

[13] Massa A, Stronati L, Aboubakr A K, et al. Numerical study of the noninertial systems: application to train coupler systems[J]. Nonlinear Dynamics, 2012, 68(1/2): 215?233.

[14] Nasr A, Mohammadi S. The effects of train brake delay time on in-train forces[J]. Proceedings of the Institution of Mechanical Engineers, Part F: Journal of Rail and Rapid Transit, 2010, 224(6): 523?534.

[15] 李亨利, 李芾, 付茂海, 等. 重載貨車坡道制動動力學及輪軌磨耗研究[J]. 鐵道科學與工程學報, 2014, 11(3): 60?64. LI Hengli, LI Fu, FU Maohai, et al. Study on the braking dynamics and wheel /rail wear for the heavy haul trucks on ramp lines[J]. Journal of Railway Science and Engineering, 2014, 11(3): 60?64.

[16] 徐倩, 王悅明, 倪純雙. 重載列車縱向沖動分布試驗研究[J]. 中國鐵道科學, 2013, 34(4): 77?83. XU Qian, WANG Yueming, NI Chunshuang. Test study on the longitudinal impulse distribution of heavy haul train[J]. China Railway Science, 2013, 34(4): 77?83.

Simulation research on the longitudinal force of heavy haul train coupler based on the improved Newmark-method

WU Chenglong, DONG Yu

(School of Automatization and Electric Engineering, Lanzhou Jiaotong University, Lanzhou 730070, China)

In the process of heavy haul train operation, the excessive longitudinal force of coupler has always been the bottleneck restricting the development of heavy haul train. The asynchronous air braking is the root of the longitudinal impulse of train, which causes the car body to squeeze the coupler to form the coupler force. The traditional method of collecting coupler force through braking characteristic test is time-consuming and labor-consuming. Therefore, in order to obtain the coupler force of heavy haul train running on different lines economically, Newmark-method was applied to the simulation analysis of coupler longitudinal force of heavy haul train. Because the train longitudinal dynamic equation was a very complex nonlinear equation, the traditional method used a lot of iterative operations to ensure the accuracy of calculation, which was time-consuming and inefficient. In this paper, the Newmark-method was improved based on the idea of increment. The nonlinear equation was dealt with directly by introducing predictive solution. Then the predictive solution was corrected. Finally the convergence approximate solution was obtained. The example shows that the improved algorithm not only ensures the calculation accuracy, but also improves the calculation efficiency. It is more suitable for the simulation calculation and analysis of the coupler longitudinal force of the long and heavy train.

heavy haul train; Newmark-method; coupler buffer; longitudinal force; prediction solution

U270.1+1

A

1672 ? 7029(2021)01 ? 0211 ? 09

10.19713/j.cnki.43?1423/u.T20200253

2020?03?31

國家自然科學基金地區科學基金資助項目(61763023)

董昱(1962?)，男，河南鄧州人，教授，從事軌道交通信息工程及控制研究；E?mail：1761853586@qq.com

(編輯 涂鵬)