交叉迭代算法求解車輛-軌道非線性耦合方程的收斂性討論

吳神花，雷曉燕

（華東交通大學鐵路環境振動與噪聲教育部工程研究中心，江西南昌330013）

交叉迭代算法求解車輛-軌道非線性耦合方程的收斂性討論

吳神花，雷曉燕

（華東交通大學鐵路環境振動與噪聲教育部工程研究中心，江西南昌330013）

摘要：運用有限元法建立車輛-軌道非線性耦合系統動力分析模型，該模型將車輛-軌道系統以輪軌接觸為界限分成車輛，軌道兩個子系統并通過輪軌接觸力的平衡和位移協調條件耦合在一起。通過交叉迭代算法分別求解車輛，軌道系統的運動方程，此時每一步都需要判斷使之滿足輪軌幾何相容條件和相互作用力平衡條件，這樣對時間步長的選取要求較高，但是如果時間步長超過某一限值，易于導致迭代失敗。引入了修正因子對輪軌接觸力進行修正，這不僅可以放寬對時間步長的選取，還能加速收斂，提高計算效率。為驗證算法的正確性，不僅進行了算例驗證，還給出了引入修正因子的交叉迭代算法求解車輛-軌道非線性耦合系統動力學方程的算例，算例中考慮了不同的時間步長和不同的修正因子對交叉迭代算法收斂速度的影響。計算結果表明引入修正因子的交叉迭代算法具有程序編制簡單、收斂速度快、用時少、精度高的優點。

關鍵詞：修正因子；時間步長；交叉迭代算法；車輛-軌道非線性耦合

作通者訊簡作介者：雷曉燕（1956—），男，教授，博士生導師，研究方向為軌道結構動力學。

隨著客運高速化、貨運重載化和運輸密度的大幅提高，車輛與軌道系統動力學問題更加突出，也日趨復雜。科研工作者對于車輛-軌道動力學理論研究也越來越深入，并取得了一系列的成果，例如翟婉明[1]建立了車輛-軌道耦合動力學統一模型；雷曉燕[2-5]采用解析的波數-頻域法建立了軌道結構單層或多層梁模型，分析了高速列車引起的軌道和大地振動；向俊,赫丹,曾慶元[6]針對無碴軌道結構提出了橫向有限條與板段單元動力分析新模型。將高速列車的動車及拖車均離散為具有二系懸掛的多剛體系統，基于彈性系統動力學總勢能不變值原理及形成系統矩陣的“對號入座”法則，建立高速列車-無碴軌道時變系統豎向振動矩陣方程，采用Wilson-θ法求解；張斌，雷曉燕等[7-8]建立列車了列車-軌道-路基耦合系統動力分析模型,提出車輛單元和軌道單元。利用車輛單元和軌道單元，考慮列車速度、路基剛度以及過渡段軌道不平順和路基剛度綜合影響因素對軌道過渡段動力特性進行分析；夏禾、曹艷梅等[9]利用解析的波數-頻域法建立了列車-軌道-大地耦合模型，將軌道-大地系統考慮為三維層狀大地上周期性支撐的歐拉梁模型，也分析了移動列車軸荷載和軌道不平順引起的動態輪軌力作用下大地的振動響應。上述關于車輛-軌道耦合動力問題的研究各有優勢，但也有其缺點：如頻域法不能解決非線性模型、非規則軌道結構模型及處理動力型不平順的問題，“對號入座法”由于車輛和軌道系統的系數矩陣將隨車輛在軌道上位置的移動而發生變化，導致在每一時步必須重新生成，使得計算量很大。而文中建立的車輛-軌道非線性耦合系統動力分析模型和引入了修正因子對輪軌接觸力進行修正的交叉迭代算法可以很好地解決這些問題。運用有限元法建立車輛-軌道非線性耦合系統動力分析模型，該模型將車輛-軌道系統以輪軌接觸為界限分成車輛，軌道兩個子系統，并通過輪軌接觸力的平衡和位移協調條件耦合在一起。通過交叉迭代算法分別求解車輛，軌道系統的運動方程，此時每一步都需要判斷使之滿足輪軌幾何相容條件和相互作用力平衡條件，這樣對時間步長的選取要求較高，但是如果時間步長超過某一限值，易于導致迭代失敗。文中引入了修正因子對輪軌接觸力進行了修正，這不僅可以放寬對時間步長的選取，還能加速收斂，提高計算效率。為證明算法的正確性，文中進行了算例驗證，還給出了引入修證因子的交叉迭代求解車輛-軌道非線性耦合系統動力學方程的算例，算例中分別考慮了不同的時間步長和不同的修正因子對交叉迭代算法收斂速度的影響。計算結果表明引入修正因子的交叉迭代算法具有程序編制簡單、收斂速度快、用時少、精度高的優點。

1　基本假設

在用有限元法建立車輛-軌道非線性耦合系統豎向振動模型時，采用以下基本假設：

1）僅考慮車輛-軌道耦合系統豎向振動效應。

2）車輛系統和軌道-路基系統沿線路方向左右對稱，可取一半結構研究。

3）上部車輛系統為附有二系彈簧阻尼的整車模型，車體和轉向架考慮沉浮振動和點頭振動。

4）輪軌間為非線性彈性接觸。

5）鋼軌被離散為二維梁單元，軌下墊板和扣件的彈性及阻尼分別用彈性系數ky1和阻尼系數cy1表示。

6）軌枕質量作為集中質量處理并僅考慮豎向振動效應；枕下道床的支承彈性系數和阻尼系數分別用ky2和cy2表示。

7）道砟質量簡化為集中質量并僅考慮豎向振動效應；道砟下路基的支承彈性系數和阻尼系數分別用ky3和cy3表示。

2　車輛子系統

車輛子系統為附有二系彈簧阻尼的整車模型，如圖1所示。圖中Mc，Jc為1/2的車體質量與轉動慣量；Mt，Jt為1/2轉向架質量與轉動慣量；ks1，ks2為車輛一、二系懸掛剛度；cs1，cs2為車輛一、二系懸掛阻尼；Mwi(i=1,2,3,4)為第i個車輪的質量；vc，c為車體沉浮振動的豎向位移、車體點頭振動的角位移；vti，ti(i=1,2)為前、后轉向架沉浮振動的豎向位移、點頭振動的角位移；vwi(i=1,2,3,4)為第i個車輪的豎向位移。考慮軌面隨機不平順，不平順幅值用表示，與四個車輪接觸處的不平順幅值分別為1，2，3，4。

圖1　上部車輛子系統模型Fig.1 Vehicle subsystem model

定義車輛單元節點位移向量：

運用Lagrange方程，可得到上部車輛子系統的動力學方程為

Pi=-Mwig+Fuli，Fuli為第i個車輪的輪軌接觸力，可根據輪軌相對接觸豎向位移由赫茲公式求得

其中：vwi，vlci分別為車輪、鋼軌在第i個輪軌接觸處的位移，i為鋼軌表面不平順，G為接觸撓度系數[10]。

3　有砟軌道子系統

有砟軌道子系統模型如圖2所示。圖3為有砟軌道單元模型，圖中v1，v4表示鋼軌的豎向位移；1，4表示鋼軌的轉角；v2，v5表示軌枕的豎向位移；v3，v6為道砟的豎向位移；mt，mb分別為四分之一的軌枕質量和兩軌枕間道砟質量；l為有砟軌道軌枕間距。

圖2　有砟軌道子系統模型Fig. 2 Ballastless track subsystem model

定義有砟軌道單元的結點位移為

圖3　有砟軌道單元模型Fig.3 Ballasted track element model

運用有限元法，可得到下部有砟軌道子系統的有限元方程為

4　引入修正因子的交叉迭代算法

采用Newmark數值積分法求解車輛-軌道非線性耦合系統動力學方程，下面給出了引入修正因子的交叉迭代算法求解車輛-軌道非線性耦合系統動力學方程的步驟。

待求解的車輛和軌道子系統方程為二階常微分方程組。

對于車輛結構，見圖1，公式(2)可改寫為

式中：kw為輪軌線性化接觸剛度；p0為車輪靜荷載。對于軌道結構，見圖2和圖3，見方程（6）

其中，ci() i=0,1,2,…,7為Newmark算法系數［11］。

現給出主要計算步驟如下：

I初始計算

II對時間步長循環

設在時間步長t，已進行了k-1次迭代，現考察第k次迭代：

定義收斂準則為

收斂準則也可以定義為力的收斂準則，但是力的收斂準則比位移收斂準則的計算速度慢。

（7）對軌道位移進行收斂性判別

（a）如果收斂性得到滿足，轉步驟II，進入下一時間步長循環，令t=t+t，并取

繼續計算，直至整個計算時域T。

（b）如果收斂性不滿足，轉步驟II，令k=k+1，進入下一迭代步循環。

5　模型驗證

為了驗證模型和方法的正確性，利用文中提出的修正后的交叉迭代算法計算文獻［7］中的算例。分別用兩種方法分析整車通過時車輛和軌道的動態響應，其中車輛為我國高速鐵道車輛，軌道為60 kg·m-1有砟軌道，列車速度為252 km·h-1，假設軌道為完全平順，其它參數見文獻［9］。計算結果見圖4，可見兩者吻合良好，說明引入修正因子的交叉迭代算法具有有效性和可行性。

圖4　鋼軌撓度曲線Fig.4 Time history of rail deflection

6　算例分析

如果不引入修正因子，時間步長很小才能收斂，而時間步長過小則會導致收斂速度慢，耗時長。

文中給出了修正因子交叉迭代求解車輛-軌道非線性耦合系統動力學方程的算例，算例中車輛參數為考慮高速整車CRH3列車，參數見表1，列車運行速度為200 km·h-1，運行時間為3 s。軌道結構為60 kg·m-1有砟軌道，無縫線路，其參數見表2。采用美國六級不平順功率譜密度函數作為軌道不平順的激勵，收斂精度均取為10e-8。由于時間步長選取大于0.002 s時，即使引入修正因子，算法仍不收斂。為了討論交叉迭代算法求解方程的收斂情況和計算效率，這里的時間步長分別取為0.002，0.001，0.000 5，0.000 25 s，修正因子則隨收斂情況適當選取，具體結果見表3-表6。表中的迭代次數約定為算完一次交叉迭代，即在某一時間步內軌道結構方程和車輛方程均完成一次求解；迭代次數按最小和最大迭代次數統計。圖5為利用力的收斂準則與利用位移收斂準則得到的輪軌力相對差值、圖6為時間步長0.000 25 s，修正因子在0.6與0.7時的輪軌力時程相對差值。

表1　CRH3高速列車車輛參數表[7]Tab.1 Vehicle parameters of CRH3

表2　我國干線60 kg/m軌道結構基本參數[8]Tab.2 Parameters for ballasted track

圖5　力收斂準則與位移收斂準則下的輪軌力時程的相對差值Fig.5 The relative difference of time history of the wheel-rail force between the force convergence and displacement convergence criterion

由圖5看出，力收斂準則與位移收斂準則下輪軌力時程的相對差值在0.6%以內，基本可以認為相等，這說明力收斂和位移收斂這兩種常用的收斂模式均可以采納，但在相同的情況下利用位移收斂準則時程序運行的時間更短，所以本文均采取位移收斂模式。

圖6　修正因子為0.6與0.7時輪軌力時程相對差值Fig.6 The relative difference of time history of the wheel-rail force between the correction factor of 0.6 and 0.7

由圖6看出，不同修正因子輪軌力相對誤差在0.3%以內，基本可以認為相等，這說明不同修正因子的引入對車輛-有砟軌道非線性耦合系統的動力響應的計算結果沒有影響。

表3　時間步長0.002 s時收斂情況和計算效率Tab.3 The convergence and computational efficiency in time step of 0.002 s

由表3可知，當時間步長為0.002 s時，修正因子大于0.1時算法不收斂，在修正因子為0.1時，程序運行時間為226 s。

表4　時間步長0.001 s時收斂情況和計算效率Tab.4 The convergence and computational efficiency in time step of 0.001 s

由表4可知，當時間步長為0.001 s時，修正因子大于0.3時算法不收斂，在程序收斂的情況下修正因子越大迭代次數越少，程序運行時間也越短，在時間步長1 ms時程序運行時間最小為199 s。

表5　時間步長0.000 5 s時收斂情況和計算效率Tab.5 The convergence and computational efficiency in time step of 0.000 5 s

由表5可知，當時間步長為0.000 5 s時，修正因子大于0.5時算法不收斂；在收斂的情況下修正因子越大迭代次數越少，程序運行時間也越短，此時程序運行時間最小為442 s。

表6　時間步長0.000 25 s時收斂情況和計算效率Tab.6 The convergence and computational efficiency in time step of 0.000 25 s

由表6可知當時間步長為0.000 25 s時，修正因子大于0.7時算法不收斂，在收斂的情況下修正因子越大迭代次數越少，程序運行時間也越少，此時程序運行時間最小為1 388 s。

7　　結論

1）若不引入修正因子，則需要選取較大的時間步長，而時間步長超過某一限值，又易導致迭代失敗。

2）若時間步長過長即使引入修正因子算法仍不收斂；在能收斂的情況下，時間步長、修正因子越大，收斂速度越快，程序運行時間越短；但時間步長較長時，取很小的修正因子時算法才能收斂。為提高算法收斂速度，需要綜合考慮時間步長和修正因子的影響。

3）將車輛-軌道非線性耦合系統分解為車輛子系統和軌道子系統，可分別獨立求解，既降低了分析問題的規模，又減小了程序設計的難度。同時，由于兩個子系統有限元方程的系數矩陣全部為定值，且為對稱矩陣，只須對方程的系數矩陣一次求逆，在每一時步的迭代中只須進行回代計算，因而極大地提高了計算速度。目前已有的基于“對號入座”求解車輛-軌道耦合系統動力響應的算法，隨著車輛在軌道上的位置不斷改變，有限元耦合方程的系數矩陣也在不斷變化，因此在每一時步的計算中，需要進行求逆運算，大大降低了計算效率。目前已有的基于“對號入座”求解車輛-軌道耦合系統動力響應的算法，隨著車輛在軌道上的位置不斷改變，有限元耦合方程的系數矩陣也在不斷變化，因此在每一時步的計算中，需要進行求逆運算，大大降低了計算效率。本文通過對列車在200 m長軌道上運行引起車輛和軌道動力響應實例的仿真分析，得到結論：采用本文提出的引入修正因子的“交叉迭代算法”在普通計算機工作站上運行用時只須199 s，而基于“對號入座”算法用時則需要60 min，一般情況下，前者的計算效率是后者的10多倍，而且求解問題的規模越大則計算效率越高。

應當指出，文中提出的算法具有普遍適用性，可廣泛應用于移動荷載作用下的各類線性和非線性問題的分析中。

參考文獻:

[1]翟婉明.車輛-軌道垂向系統的統一模型及其耦合動力學原理[J].鐵道學報, 1992,14(3):10-21.

[2]雷曉燕.軌道臨界速度與軌道強振動研究[J].巖土工程學報,2006,28(3):419-422.

[3]雷曉燕.軌道結構動力分析的傅里葉變換法[J].鐵道學報,2007,29(3):67-71.

[4]雷曉燕.高速列車誘發地面波與軌道強振動研究[J].鐵道學報,2006,28(3):78-82.

[5]雷曉燕.高速鐵路軌道振動與軌道臨界速度的傅里葉變換法[J].中國鐵道科學,2007,28(6):30-34.

[6]向俊,赫丹,曾慶元.橫向有限條與無砟軌道板段單元的車軌系統豎向振動分析法[J].鐵道學報, 2007,29(4):64-69

[7]雷曉燕,張斌,劉慶杰.軌道過渡段動力特性的有限元分析[J].中國鐵道科學,2009,30(5):15-21.

[8]雷曉燕,張斌,劉慶杰.列車-軌道系統豎向動力分析的車輛軌道單元模型[J].振動與沖擊,2010,29(3):168-173.

[9] H XIA, Y M CAO. Theoretical modeling and characteristic analysis of moving-train induced ground vibrations [J]. Journal of Sound and Vibration,2010,329(7):819-5163.

[10]雷曉燕.軌道力學與工程新方法[M].北京:中國鐵道出版社,2000:45-46.

[11]雷曉燕.有限元法[M].北京:中國鐵道出版社,2000:177-180.

（責任編輯王建華）

Convergence Condition of Cross Iterative Algorithm for Vehicle-Track Nonlinear Coupling Equations

Wu Shenhua, Lei Xiaoyan
(Engineering Research Center of Railway Environment Vibration and Noise of the Ministry of Education, East China Jiaotong University, Nanchang 330013, China)

Abstract:A model for dynamic analysis of vehicle-track nonlinear coupling system was established by finite ele?ment method. The whole system was divided into two subsystems, i.e., the vehicle subsystem considered as second?ary spring vehicle model and the track subsystem regarded as three discrete elastic beams model. Coupling of the two systems was achieved by equilibrium conditions for wheel-rail nonlinear contact forces and geometrical com?patibility conditions. A cross-iterative algorithm was presented to solve the dynamics equations of vehicle-track nonlinear coupling system, but every step must be judged to satisfy equilibrium conditions for wheel-rail nonlinear contact forces and geometrical compatibility conditions, leading to the selection of a higher time step that if the time step exceeded a certain value it would be easy to cause failure of the iteration. In view of this situation, a re?laxation technique was introduced to correct the wheel-rail contact force, which not only could broaden the selec?tion of the time step but also could accelerate the iterative convergence rate. By contrast with the example of refer?ence, the correctness of the algorithm was verified. Examples of cross-iterative algorithm of the relaxation tech?nique were given, in which the influence of different time steps and relaxation techniques were considered. Results demonstrated that the cross-iterative algorithm had the advantages of simple programming, fast convergence rate, less computation time and high accuracy.

Key words:relaxation technique; time step; cross-iterative algorithm; vehicles-track nonlinear coupling system

作者簡介：吳神花（1989—），女，碩士研究生，研究方向為軌道結構動力學。

基金項目：國家自然科學基金項目(U1134107)

收稿日期：2014-12-25

文章編號：1005-0523（2015）03-0023-09

中圖分類號：U213.2+12

文獻標志碼：A