秦佳良，劉林芽，2

(1.華東交通大學(xué) 軌道交通基礎(chǔ)設(shè)施性能監(jiān)測與保障國家重點實驗室，江西 南昌 330013；2.萍鄉(xiāng)學(xué)院 工程與管理學(xué)院，江西 萍鄉(xiāng) 337055)

鐵路交通是國家發(fā)展的重大公共基礎(chǔ)設(shè)施,是現(xiàn)代化綜合交通運輸體系的基礎(chǔ)骨干,對促進經(jīng)濟社會的快速發(fā)展起到了十分顯著的作用。我國在進入21世紀以來,在高速鐵路、重載鐵路、城市軌道交通等領(lǐng)域取得了世界矚目的成就。軌道交通通過輪軌相互作用來實現(xiàn)其功能,保持軌道良好的平順性,降低機車車輛對軌道和軌下基礎(chǔ)的動力作用,減少養(yǎng)護維修的成本,以及新型機車車輛和新型軌道結(jié)構(gòu)的設(shè)計、制造等問題的解決都需要對車輛與軌道系統(tǒng)之間相互作用的動力特性展開深入研究[1]。因此準確、高效地求解車輛-軌道耦合系統(tǒng)的動態(tài)響應(yīng)對保證軌道交通的可持續(xù)發(fā)展具有重大的理論和現(xiàn)實意義。

車輛-軌道耦合系統(tǒng)的動力響應(yīng)模型求解一直是車輛動力學(xué)和軌道動力學(xué)研究中最關(guān)鍵的基礎(chǔ)問題[2]。近幾十年以來,國內(nèi)外眾多學(xué)者曾先后對車輛-軌道耦合作用引起的機車車輛、軌道結(jié)構(gòu)、路基、隧道、橋梁及大地振動等方向展開了深入的研究[3-7],提出一系列繁簡各異的輪軌動力學(xué)理論、算法、模型、程序及軟件,部分實現(xiàn)了工程應(yīng)用,有效地解決了在軌道交通可持續(xù)發(fā)展中面臨的一些瓶頸問題,也取得了一些較為顯著的研究成果[8-12]。在一般應(yīng)用情況下,車輛-軌道耦合系統(tǒng)的動力學(xué)分析模型較復(fù)雜且動力自由度數(shù)目也較多,采用解析法直接分析和求解其動力學(xué)微分方程組會十分復(fù)雜且困難,因此目前一般都普遍采用直接數(shù)值積分法來求解。直接數(shù)值積分法分為顯式數(shù)值積分法和隱式數(shù)值積分法兩種不同的形式。隱式數(shù)值積分法一般具有穩(wěn)定性較好、計算的精度較高等特點,普遍使用的隱式數(shù)值積分法一般包括Newmark-β法、Wilson-θ法、Park法和Hobolt法等。但隱式數(shù)值積分法每步積分均需求解大型代數(shù)方程組,這對計算規(guī)模大的工程問題帶來了巨大的挑戰(zhàn)。與隱式數(shù)值積分法不同,顯式數(shù)值積分法則一般具有計算過程較簡單、計算效率較高等一些優(yōu)點,普遍使用的顯式數(shù)值積分法一般包括新型快速顯式積分法、四階Runge-Kutta法和中心差分法等。顯式數(shù)值積分法明顯的缺點就是其穩(wěn)定性或計算精度一般比隱式數(shù)值積分法略差?？梢?顯式數(shù)值積分與隱式數(shù)值積分的優(yōu)缺點是互補的。鑒于此,我們在求解車輛-軌道非線性耦合系統(tǒng)振動方程時,嘗試將顯式數(shù)值積分與隱式數(shù)值積分聯(lián)合使用,以期獲得比較好的綜合性能。

本文利用有限元法分別建立了車輛子系統(tǒng)和軌道子系統(tǒng)的動力學(xué)分析模型,然后采用顯式數(shù)值積分法和隱式數(shù)值積分法分別求解車輛子系統(tǒng)和軌道子系統(tǒng)的動力學(xué)方程,提出了車輛-軌道非線性耦合系統(tǒng)振動分析的分離迭代和分離同步兩種求解方法,并對兩種方法進行了算例驗證。最后通過進一步研究時間步長對輪軌耦合系統(tǒng)振動響應(yīng)的影響,以及對比不同時間步長下車輛-軌道非線性耦合系統(tǒng)的數(shù)值計算效率,分析兩種算法在計算穩(wěn)定性、準確性和數(shù)值編程求解效率等方面的差異,為求解車輛-軌道非線性耦合相互作用問題提供可供選擇的高效數(shù)值計算方法。

1 車輛-軌道非線性耦合系統(tǒng)動力學(xué)方程

本文建立車輛-軌道非線性耦合系統(tǒng)分析模型時,由于車輛子系統(tǒng)和軌道子系統(tǒng)沿線路方向左右對稱,因此本文只取半結(jié)構(gòu)展開研究。且本文僅考慮輪軌豎向動力效應(yīng),輪軌間考慮為Hertz非線性接觸。

1.1 車輛子系統(tǒng)動力學(xué)方程

車輛子系統(tǒng)可以考慮為一個10自由度的具有一、二系懸掛的整車模型,模型分別考慮車體、轉(zhuǎn)向架的沉浮和點頭運動以及輪對的沉浮運動,具體情況見圖1。圖1中,Mc、Jc分別為車體的質(zhì)量、轉(zhuǎn)動慣量;Mt、Jt分別為轉(zhuǎn)向架的質(zhì)量、轉(zhuǎn)動慣量;Mwi(i=1,2)為第i個車輪的質(zhì)量;Ks1、Cs1分別為車輛的一系懸掛剛度、阻尼;Ks2、Cs2分別為車輛的二系懸掛剛度、阻尼;νc、νti(i=1,2)分別為車體、前后轉(zhuǎn)向架沉浮運動的豎向位移;θc、θti(i=1,2)分別為車體、前后轉(zhuǎn)向架點頭運動的角位移;νwi(i=1,2,3,4)為第i個車輪沉浮運動的豎向位移;l1、l2分別為車輛固定軸距之半、車輛定距之半。

車輛子系統(tǒng)的位移向量XV為

XV={vcθcvt1θt1vt2θt2vw1vw2vw3vw4}T

(1)

由Lagrange方程我們可以直接得到車輛子系統(tǒng)的動力學(xué)方程為

(2)

QV=-g{Mc0Mt0Mt0MwMwMwMw}T

(3)

FVT={0 0 0 0 0 0F1F2F3F4}T

(4)

式中:g為重力加速度;Fi(i=1,2,3,4)為輪軌相互用力,可用Hertz非線性彈性接觸理論來計算,其表達式為

(5)

其中,νlci、ηi(i=1,2,3,4)分別為第i個輪軌接觸處的鋼軌位移、軌道的不平順幅值;G為輪軌接觸常數(shù)[5]。

1.2 軌道子系統(tǒng)動力學(xué)方程

以路基上采用的CRTS Ⅱ型板式無砟軌道結(jié)構(gòu)為例,建立模擬板式無砟軌道子系統(tǒng)的動力學(xué)分析模型。針對CRTS Ⅱ型板式無砟軌道子系統(tǒng)結(jié)構(gòu)形式的特點,將軌道子系統(tǒng)沿線路縱向取半結(jié)構(gòu)進行分析研究,并利用有限元法來模擬板式無砟軌道子系統(tǒng)。將沿線路縱向兩相鄰扣件之間的軌道結(jié)構(gòu)考慮為一個板式無砟軌道單元,單元長度記為l,單元模型見圖2。模型中,將鋼軌模擬為離散黏彈性點支承梁單元,Ky1、Cy1分別為扣件的彈性系數(shù)、阻尼系數(shù);將軌道板和底座板模擬為連續(xù)黏彈性支承梁單元,Ky2、Cy2分別為CA砂漿層的彈性系數(shù)、阻尼系數(shù),Ky3、Cy3分別為路基層的彈性系數(shù)、阻尼系數(shù)。圖2中,ν1、ν4分別為鋼軌兩端的豎向位移;θ1、θ4分別為鋼軌兩端的轉(zhuǎn)角;ν2、ν5分別為軌道板兩端的豎向位移;θ2、θ5分別為軌道板兩端的轉(zhuǎn)角;ν3、ν6分別為底座板兩端的豎向位移;θ3、θ6分別為底座板兩端的轉(zhuǎn)角。

圖2 板式無砟軌道單元模型

(6)

(7)

最后利用有限元法的“對號入座”法則,根據(jù)板式無砟軌道單元的質(zhì)量矩陣、阻尼矩陣和剛度矩陣,組集整個軌道子系統(tǒng)的質(zhì)量矩陣、阻尼矩陣、剛度矩陣,可得到軌道子系統(tǒng)的動力學(xué)方程為

(8)

2 顯隱式求解算法計算流程

在用數(shù)值方法求解車輛-軌道非線性耦合系統(tǒng)分析模型的動力學(xué)方程時,因車輛子系統(tǒng)的質(zhì)量矩陣為對角陣,數(shù)值計算時可以采用新型顯式積分法來求解車輛子系統(tǒng)的動力學(xué)方程;軌道子系統(tǒng)采用有限元法建立,其質(zhì)量矩陣為非對角陣,數(shù)值計算時可以采用Newmark-β積分法來求解軌道子系統(tǒng)的動力學(xué)方程。分別建立車輛-軌道非線性耦合系統(tǒng)振動分析的分離迭代和分離同步兩種求解方法。分離迭代法是將車輛-軌道非線性耦合系統(tǒng)劃分為車輛和軌道兩個非時變的子系統(tǒng)模型,通過輪軌位移協(xié)調(diào)或輪軌力平衡關(guān)系將兩個子系統(tǒng)模型聯(lián)系起來,兩個子系統(tǒng)模型之間通過迭代計算來分別求解車輛子系統(tǒng)模型和軌道子系統(tǒng)模型的動力學(xué)方程。分離同步法其實也是將車輛-軌道非線性耦合系統(tǒng)分析模型劃分為車輛和軌道兩個非時變子系統(tǒng)模型,但車輛子系統(tǒng)模型和軌道子系統(tǒng)模型的動力學(xué)方程是同步求解的,不需要進行迭代求解。

2.1 顯隱式積分求解格式

(9)

(10)

式中:Δt為時間積分步長;ψ、φ分別為控制積分方法特性的兩個獨立參數(shù),計算起步時只需令ψ=φ=0。

(11)

(12)

(13)

式中:c0=1/(αΔt2);c1=δ/(αΔt);c2=1/(αΔt);c3=1/(2α)-1;c4=δ/α-1;c5=Δt(δ/α-2)/2;c6=Δt(1-δ);c7=δΔt;Newmark-β積分常數(shù)α和δ分別取0.25和0.5。

2.2 顯隱式分離迭代法求解流程

在開始求解時先假設(shè)車輛子系統(tǒng)和軌道子系統(tǒng)從靜止起步,即初始時刻t=0時,車輛子系統(tǒng)和軌道子系統(tǒng)的位移響應(yīng)、速度響應(yīng)和加速度響應(yīng)均為零。在對時間步長進行循環(huán)計算時,假設(shè)在t+Δt時刻,兩個子系統(tǒng)動力學(xué)方程求解已進行k次迭代,現(xiàn)在考察第k+1次迭代。

當(dāng)k為偶數(shù)時,令

(14)

并計算差值Δ1(i)為

(15)

當(dāng)k為奇數(shù)時,令

(17)

注意到Aitken加速法是每迭代兩次才改善一次。

Step6根據(jù)軌道子系統(tǒng)位移差值Δ進行收斂性判別,即

(18)

收斂條件為

(19)

Step7如果計算滿足收斂性,則可以轉(zhuǎn)入下一時間步長后繼續(xù)循環(huán)計算,直至整個時域T。如果計算不滿足收斂性,則令k=k+1,進入下一迭代步后繼續(xù)循環(huán)計算。

為進一步提高顯隱式分離迭代法的計算效率,參照文獻[14]中的改進方法,根據(jù)式( 9 )、式(10)和式(11),可以在迭代步循環(huán)計算前,先計算車輛子系統(tǒng)的位移和速度,以及軌道子系統(tǒng)位移響應(yīng)在迭代步循環(huán)過程中的不變量,在迭代求解過程中直接調(diào)用即可,這樣可以減少計算工作量。

2.3 顯隱式分離同步法求解流程

與顯隱式分離迭代法一樣,顯隱式分離同步法在開始求解時車輛子系統(tǒng)和軌道子系統(tǒng)也從靜止起步,即假設(shè)初始時刻t=0時,車輛子系統(tǒng)和軌道子系統(tǒng)的位移響應(yīng)、速度響應(yīng)和加速度響應(yīng)均為零。現(xiàn)假設(shè)已經(jīng)求得t-Δt和t時刻車輛子系統(tǒng)和軌道子系統(tǒng)的位移響應(yīng)、速度響應(yīng)和加速度響應(yīng),現(xiàn)在需要求解t+Δt時刻兩個子系統(tǒng)的位移響應(yīng)、速度響應(yīng)和加速度響應(yīng)。主要計算步驟如下:

3 顯隱式求解算法驗證

為驗證本文提出的顯隱式分離迭代法和顯隱式分離同步法的正確性,與文獻[15]中計算的算例進行比較。計算條件為假設(shè)列車編組為1動+1拖的高速列車在博格板式無砟軌道上運行,運行速度為200 km/h,軌道高低不平順激振源考慮成波長為20 m、波幅為6 mm的周期性正弦函數(shù),其計算結(jié)果見圖3,文獻[15]的計算結(jié)果見圖4。

圖3 本文計算結(jié)果

圖4 參考文獻[15]中計算結(jié)果

由圖3和圖4的計算結(jié)果對比可知,鋼軌和博格板的位移響應(yīng)波形均符合物理概念,與文獻[15]計算得到的響應(yīng)幅值與變化規(guī)律基本一致,證明了本文提出的顯隱式分離迭代法和顯隱式分離同步法的正確性和可行性。

4 顯隱式求解算法收斂性對比分析

本文為對比分析顯隱式分離迭代法和顯隱式分離同步法的收斂性,以CRH3型高速動車通過CRTSⅡ型板式無砟軌道為計算算例。其中車輛的運行速度設(shè)置為300 km/h,板式無砟軌道單元的長度取0.65 m,在進行數(shù)值計算時軌道子系統(tǒng)的長度設(shè)為360個軌道單元的長度。軌道高低不平順設(shè)置為車輛-軌道非線性耦合系統(tǒng)的激勵源,將其考慮成波長為12.5 m、波幅為3 mm的周期性正弦函數(shù)。車輛子系統(tǒng)和軌道子系統(tǒng)的結(jié)構(gòu)參數(shù)按照文獻[16]中的參數(shù)取值。

4.1 時間步長的確定

在應(yīng)用中最為重要的問題就是時間步長的確定,這直接影響到數(shù)值計算是否收斂以及計算效率。我們采用數(shù)值試驗方法,對各種時間步長(最小時間步長為10-3ms)進行了數(shù)值模擬,計算結(jié)果見圖5、圖6。

圖5 輪軌振動加速度的數(shù)值模擬結(jié)果

圖6 輪軌力的數(shù)值模擬結(jié)果

由圖5和圖6可知,時間步長對輪對加速度和輪軌力的影響相對較小,而對鋼軌加速度的影響較大。顯隱式分離迭代法的臨界收斂時間步長為1 ms,超過1 ms后計算結(jié)果將不收斂,最大有效時間步長可取為0.5 ms。顯隱式分離同步法的臨界收斂時間步長和最大有效時間步長均為0.1 ms,超過0.1 ms后計算結(jié)果將不收斂。

4.2 計算效率的對比分析

為對比顯隱式分離迭代法和顯隱式分離同步法在計算效率上的差異,在計算機上運行程序,計算時間步長分別取0.5、0.1、0.05、0.01 ms時,考察兩種方法的計算耗時。不同時間步長下車輛-軌道非線性耦合系統(tǒng)的求解時間見表1。

表1 不同時間步長耦合系統(tǒng)求解時間

由表1可知,采用顯隱式分離迭代法求解時,當(dāng)時間步長取得較大時,采用Aitken加速法可以增強系統(tǒng)迭代求解的計算穩(wěn)定性,但當(dāng)時間步長取得越小,Aitken加速法的作用會越來越小。這是因為當(dāng)時間步長較大時,耦合系統(tǒng)在下一時刻的平衡狀態(tài)會和上一時刻耦合系統(tǒng)的平衡狀態(tài)產(chǎn)生嚴重偏離的情況,這將導(dǎo)致耦合系統(tǒng)數(shù)值解的發(fā)散而使得計算結(jié)果無法收斂,或者耦合系統(tǒng)需要多次迭代才能使計算結(jié)果收斂,此時采用Aitken加速法可以減少迭代次數(shù),可以起到增強迭代計算穩(wěn)定性的作用。而當(dāng)時間步長取得越小時,耦合系統(tǒng)在下一時刻的平衡狀態(tài)會和上一時刻耦合系統(tǒng)的平衡狀態(tài)越接近,這會使得耦合系統(tǒng)非常容易地再次達到其平衡位置,因此耦合系統(tǒng)的收斂速度會較快且Aitken加速法的作用會越來越小。

當(dāng)采用相同時間步長時,顯隱式分離同步法的計算效率要比顯隱式分離迭代法更高。這是因為采用顯隱式分離迭代法求解時,每一時間步內(nèi)車輛子系統(tǒng)和軌道子系統(tǒng)的動力學(xué)方程至少分別求解2次才能收斂,而采用顯隱式分離同步法求解時每一時間步內(nèi)只需求解一次耦合系統(tǒng)的動力學(xué)方程,因此顯隱式分離同步法的計算效率會更高。但是顯隱式分離迭代法可以通過選取較大的時間步長來提高耦合系統(tǒng)的計算效率,在實際工程應(yīng)用時可根據(jù)實際情況選用合適的求解方法。

5 結(jié)論

1)提出了基于顯隱式積分格式的車輛-軌道非線性耦合系統(tǒng)分離迭代和分離同步兩種數(shù)值計算方法,并通過算例驗證了兩種算法的正確性。

2)對車輛-軌道非線性耦合系統(tǒng)進行數(shù)值分析時,顯隱式分離迭代法和顯隱式分離同步法的最大有效時間步長分別是0.5、0.1 ms。

3)當(dāng)時間步長較大時,采用Aitken加速法可以起到增強顯隱式分離迭代法的計算穩(wěn)定性的作用。但隨著時間步長的減小,Aitken加速法的作用會越來越小。

4)相同時間步長下顯隱式分離同步法的計算效率要比顯隱式分離迭代法高,但是可以采用較大的時間步長來提高顯隱式分離迭代法的計算效率。