車輛-軌道-橋梁豎直耦合振動程序設(shè)計及仿真

雷震宇 閆 旭

(同濟大學(xué)鐵道與城市軌道交通研究院，201804，上?！蔚谝蛔髡?，副教授)

軌道交通車輛在高架橋上運行時產(chǎn)生的結(jié)構(gòu)振動和噪聲問題一直困擾著人們的生活。輪軌復(fù)雜的非線性接觸是該問題的核心。在輪軌關(guān)系的仿真分析中，應(yīng)用最廣泛的兩款軟件為 ADAMS 和SIMPACK AG。雖然，兩者的輪軌模塊可以很好地體現(xiàn)車輛和軌道的耦合作用，但對需要考慮柔性基礎(chǔ)(橋梁)的剛?cè)峄旌夏Ｐ偷姆抡?，兩款軟件都需要通過與外部有限元軟件的接口實現(xiàn)對復(fù)雜柔性體的支持，不但操作繁瑣，而且參數(shù)設(shè)置不當(dāng)還可能會影響柔性體傳輸?shù)恼_性，使得仿真結(jié)果的可參考性很差［1］。鑒于現(xiàn)有軟件對車輛－軌道－橋梁耦合系統(tǒng)的仿真功能還有待完善，許多學(xué)者都通過直接編程的方式對車橋耦合系統(tǒng)的動力特性進行了相關(guān)研究［2－9］。但是，前人的研究多數(shù)沒有考慮軌道的參振作用，或是車輛和軌道的模型過于簡化而無法準(zhǔn)確得到軌道結(jié)構(gòu)的振動特性。本文從耦合系統(tǒng)的振動微分方程出發(fā)，建立了完整的車－線－橋豎向耦合模型，并以工程中常用的橋上彈性支承塊式無砟軌道為例，利用MATLAB 軟件編寫了求解耦合系統(tǒng)豎直振動時變特性的程序;考慮到系統(tǒng)各組成部件動力特性的差異和部件之間耦合作用的強弱，對車輛、軌道和橋梁分別予以合理簡化，在程序編寫時將系統(tǒng)化整為零，將大型非線性時變微分方程組分塊求解，降低了整體求解對計算機內(nèi)存和速度的要求，節(jié)省了計算成本。仿真結(jié)果顯示，簡化模型和計算程序可以保證分析的準(zhǔn)確性。

1 豎向耦合振動模型的建立

1.1 物理模型的簡化

對物理模型作了如下的簡化:

1)車輛采用具有二系懸掛的客車模型［10］。一系和二系懸掛簡化為線性彈簧和粘性阻尼;輪對、構(gòu)架和車體簡化為7 個剛體，考慮輪對的沉浮運動、構(gòu)架和車體的沉浮和點頭運動，共有10 個自由度;輪對與鋼軌的接觸用Hertz 非線性彈簧模擬。

2)軌道選用彈性支承塊式無砟軌道。鋼軌視為連續(xù)彈性離散點支撐的Euler 梁，只考慮豎向彎曲振動;扣件和膠墊的作用通過線性彈簧和粘性阻尼模擬;將彈性支承塊視為剛性質(zhì)量塊，只考慮豎向振動;橋梁跨度范圍之外的支承塊用彈簧和阻尼與剛性基礎(chǔ)連接;將道床板和混凝土底座視為橋梁的一部分，考慮軌道高低不平順的影響。

3)橋梁簡化成線性等截面簡支梁，只考慮平面內(nèi)的彎曲振動，道床塊與橋面的相互作用通過線性彈簧和粘性阻尼模擬，忽略橋梁支座和橋墩的變形。

經(jīng)過簡化后的車輛－軌道－橋梁豎向振動系統(tǒng)如圖1所示。

圖1 車輛－軌道－橋梁豎向耦合振動簡化模型

1.2 耦合振動方程的推導(dǎo)

耦合系統(tǒng)振動方程的建立方法主要有分離子系統(tǒng)法和整體法兩種。兩種方法各有其優(yōu)缺點。分離子系統(tǒng)法能夠體現(xiàn)子系統(tǒng)之間復(fù)雜的相互作用力和位移協(xié)調(diào)關(guān)系，可以對子系統(tǒng)單獨求解，且系統(tǒng)矩陣不具時變性，減小了計算工作量，但建模較為繁瑣［11］;整體法建模思路清晰，可以避開輪軌關(guān)系的不確定性，但需要求解大型時變非線性方程組而非常耗時。為研究輪軌之間的相互作用，同時為求解方便，本文利用子系統(tǒng)建模法建立系統(tǒng)振動方程，以輪軌接觸面和道橋接觸面為界，將耦合系統(tǒng)分為車輛子系統(tǒng)、軌道子系統(tǒng)和橋梁子系統(tǒng)，各子系統(tǒng)之間通過位移相容條件和靜力平衡條件進行耦合。其方法如下。

1.2.1 車輛子系統(tǒng)的豎向振動方程

根據(jù)達(dá)朗貝爾原理，可得多剛體組成的車輛系統(tǒng)的動平衡方程［10］，寫成矩陣的表達(dá)形式為

式中:

Mv——車輛子系統(tǒng)的質(zhì)量矩陣;

Kv——車輛子系統(tǒng)的剛度矩陣;

Cv——車輛子系統(tǒng)的阻尼矩陣;

Fr－v——車輛子系統(tǒng)的廣義荷載向量;

uv——車輛子系統(tǒng)的廣義位移向量。

1.2.2 軌道子系統(tǒng)的豎向振動方程

鋼軌的振動方程為4 階偏微分方程，采用振型疊加法［12］得到正則化坐標(biāo)下的振動方程為

式中:

qk——鋼軌 k 階模態(tài)坐標(biāo);

lr——鋼軌計算長度;

pwrj——第 j號輪對的輪軌力;

Ns——計算范圍內(nèi)鋼軌支點總數(shù);

Zk(x)——鋼軌振型函數(shù);

xsi——第i 個鋼軌支點的位置;

xwj——第 j號輪對的位置。

根據(jù)不同位置處彈性支承等效質(zhì)量塊的受力狀態(tài)分別建立Ns個振動方程，并將鋼軌和等效質(zhì)量塊的振動方程在廣義坐標(biāo)下進行耦合，最終得到矩陣表達(dá)形式的軌道子系統(tǒng)的振動方程為

式中:

Mrs，Krs，Crs，urs——分別為軌道子系統(tǒng)的質(zhì)量矩陣、剛度矩陣、阻尼矩陣和廣義位移向量;

Fv－rs，F(xiàn)br－rs——分別為車輛和橋梁對軌道的廣義荷載向量。

1.2.3 橋梁子系統(tǒng)的豎向振動方程

用平面梁單元將橋梁沿跨度均分，保證在彈性支承塊底部都有對應(yīng)橋面節(jié)點。選定梁單元的位移模式為線性函數(shù)，單元質(zhì)量矩陣采用一致質(zhì)量矩陣，利用最小勢能原理和坐標(biāo)轉(zhuǎn)換矩陣得到整體坐標(biāo)下梁單元的單元剛度和質(zhì)量矩陣;將單元剛度矩陣和質(zhì)量矩陣進行擴展，集合成整體剛度和質(zhì)量矩陣。整體阻尼矩陣采用Rayleigh 阻尼，加入節(jié)點力后可得到橋梁豎向振動的有限元矩陣方程

式中，Mbr，Kbr，Cbr，F(xiàn)rs－br，ubr分別為橋梁子系統(tǒng)的質(zhì)量矩陣、剛度矩陣、阻尼矩陣、廣義荷載向量和廣義位移向量。

1.3 耦合方程的數(shù)值解法

3 個微分方程組式(1)、式(3)和式(4)的右端荷載項都包含相鄰子系統(tǒng)的位移向量，體現(xiàn)了方程之間的耦合關(guān)系。借助數(shù)值方法在時域內(nèi)對3 個子系統(tǒng)分別進行求解。采用顯－隱式相結(jié)合的方法，在每一個時間步內(nèi)，采用新型快速顯式積分法［10］求解車輛子系統(tǒng)和軌道子系統(tǒng);采用隱式積分法中的Newmark－β 積分法［13］求解橋梁子系統(tǒng)，最終可以得到系統(tǒng)的時域振動響應(yīng)。

2 分離子系統(tǒng)法程序設(shè)計流程

直接積分法中的顯式積分法是條件收斂的，時間步長和時間步的選取對程序的穩(wěn)定性和計算精度影響很大［14］。一般要求的積分時間步長要小于所感興趣的系統(tǒng)最小振動周期的1/10［15］，并且輪軌之間的非線性接觸導(dǎo)致耦合方程的求解對初值的依賴性很大，呈現(xiàn)出“混沌”特性。在輪軌脫離的瞬間，鋼軌的模態(tài)加速度很大，對模態(tài)響應(yīng)進行顯式積分后的誤差累積很可能導(dǎo)致計算發(fā)散。因此，在程序設(shè)計時，通過對鋼軌與支承塊接觸點的位移和速度進行控制，反推鋼軌模態(tài)坐標(biāo)，可以保證積分的穩(wěn)定性。經(jīng)過調(diào)試，綜合考慮到計算效率，最終決定的積分時間步長為0.000 2 s，可以精確分析系統(tǒng)500 Hz 以內(nèi)的振動情況。分離子系統(tǒng)法程序設(shè)計流程如圖2所示。

圖2 分離子系統(tǒng)法程序設(shè)計流程

3 程序校核算例

3.1 模型參數(shù)

車輛選取 HSC 型車輛，具體參數(shù)詳見文獻(xiàn)［10］，軌道和橋梁參數(shù)取值如表1所示。

表1 軌道和橋梁結(jié)構(gòu)參數(shù)

由于橋梁單元大小對耦合系統(tǒng)動態(tài)響應(yīng)的影響很?。?6］，橋梁劃分為100 個單元進行計算。選取美國六級譜軌道高低不平順為系統(tǒng)激勵源［17］，利用程序繪制軌道不平順樣本曲線，以分析橋梁振動特性，描繪列車從進入橋梁到駛出橋梁過程中的鋼軌、支承塊及橋梁在不同時刻的變形狀態(tài)和系統(tǒng)關(guān)鍵點的位移響應(yīng)時程曲線，驗證程序的可靠性。

3.2 軌道不平順模擬及精度驗證

采用三角級數(shù)法生成軌道不平順樣本，采樣間距為0.1 m。鋼軌計算長度內(nèi)共采集樣本1 200 個，取空間頻率范圍為 0.01 ～0.5 Hz，頻率帶寬為 0.01 Hz，則生成的軌道豎向不平順樣本如圖3所示。

圖3 美國六級譜軌道高低不平順模擬曲線

利用平均周期圖法對模擬的不平順數(shù)據(jù)進行離散傅里葉變換，求得不平順樣本的功率譜密度，并與理論譜密度曲線進行對比。圖4 顯示了模擬曲線的譜密度與理論要求的譜密度在目標(biāo)空間頻率內(nèi)基本吻合，說明模擬的有效性。

3.3 橋梁振動特性

利用MATLAB 軟件編寫了有限元法求解簡支橋梁前6 階豎向無阻尼自由振動頻率和振型的程序［18］，結(jié)果如表 2所示。由表 2 可以看出，橋梁主要振動頻率在200 Hz 以內(nèi)。

圖4 美國六級譜軌道高低不平順功率譜密度曲線

3.4 程序的有效性驗證

以鋼軌左端為整體坐標(biāo)原點，行車方向為x 軸，坐標(biāo)原點x=0，橋梁處在x=45 ～75 m 范圍內(nèi)。初始時刻列車、軌道和橋梁都處于各自的靜平衡位置，且處于靜止?fàn)顟B(tài);列車前轉(zhuǎn)向架前輪處于x =30 m 處，輪軌相互接觸但無擠壓。不計初始變形，給列車水平向初速度作為系統(tǒng)初始激勵，所有仿真結(jié)果都基于此初始條件。仿真距離取80 m，自x=40 m 處開始考慮軌道不平順，為防止軌道高差突變對程序收斂性的影響，在 x=35 ～40 m 范圍內(nèi)引入 5 m 的緩和曲線，以保證軌道不平順線性緩慢增加。

假定車輛以50 m/s 的速度運行，前轉(zhuǎn)向架前輪在 0.3 s、0.6 s 和 0.9 s 到達(dá)橋頭、跨中和橋尾，后轉(zhuǎn)向架后輪在 0.71 s、1.01 s 和 1.31 s 到達(dá)橋頭、跨中和橋尾，選取這6 個時刻的狀態(tài)進行參考，按照時間先后順序繪制鋼軌、支承塊和橋梁的位移狀態(tài)如圖5a)～5f)所示。圖像顯示，程序可以有效描繪耦合系統(tǒng)組件的狀態(tài)變化。

當(dāng)車輛以50 m/s 的速度經(jīng)過高架橋區(qū)段時，鋼軌中點、中間支承塊和橋梁跨中點的時程和頻譜(對數(shù)坐標(biāo))曲線如圖6所示。

由圖6a)可知，鋼軌和支承塊的時程曲線各有4個峰值，分別對應(yīng)車輛的4 個輪對剛好經(jīng)過該點的時刻;車輪作用點處的鋼軌位移很大，作用點附近的鋼軌位移很小，甚至為0，響應(yīng)波形符合客觀實際;橋梁時程曲線只有2 個峰值，這是由橋梁振動的低頻特性造成的，在轉(zhuǎn)向架的前后輪先后經(jīng)過橋梁跨中點的時間段內(nèi)，橋梁撓度來不及回彈，2 個峰值分別對應(yīng)車輛的前后轉(zhuǎn)向架剛好經(jīng)過橋梁跨中點的時刻。

由圖6b)可知，橋梁振動頻率分布在0 ～200 Hz，主要頻率成分集中在5 Hz 以內(nèi)。這說明橋梁以一階彎曲振動為主，橋梁實際振動頻率比無阻尼振動頻率低。而鋼軌和支承塊的振動頻率分布范圍較寬，在0～500 Hz，且以低于20 Hz 的低頻振動為主;鋼軌和支承塊的頻率分布大致相同，100 ～200 Hz 的頻率成分所占比重比橋梁在同頻段內(nèi)的比重大。

圖5 鋼軌、支承塊和橋梁不同時刻的位移狀態(tài)圖像

4 仿真結(jié)果探討

4.1 軌道不平順對系統(tǒng)動力響應(yīng)的影響

車輛速度取為50 m/s，當(dāng)軌道完全平順和軌道不平順時，輪軌相互作用力、鋼軌中點豎向加速度和前轉(zhuǎn)向架前輪的豎向加速度時程曲線如圖7 和圖8所示。圖7、圖8 中時間軸從0.2 s 開始顯示，原因是程序在開始運行階段，對初始條件需要經(jīng)過多次迭代以達(dá)到系統(tǒng)穩(wěn)定，因此該時間段內(nèi)的數(shù)據(jù)沒有實際參考意義，故此舍去。

圖6 鋼軌中點、中間支承塊和橋梁跨中點的時頻圖

圖7 平順軌道的系統(tǒng)動力響應(yīng)

圖8 不平順軌道的系統(tǒng)動力響應(yīng)

圖7 和圖8 的對比說明，軌道不平順對系統(tǒng)動態(tài)特性影響很大。軌道平順性良好時，輪軌作用力在軸重附近有規(guī)律地上下波動，其最大值為156.804 kN，鋼軌中點和前轉(zhuǎn)向架前輪的豎向加速度峰值都很小。此時車輛、軌道和橋梁的動力作用不明顯，輪軌密貼接觸，列車運行十分平穩(wěn)。軌道不平順時，輪軌作用力仍圍繞軸重上下波動，但波動劇烈、雜亂無章，多次出現(xiàn)輪軌力為0，即某車輪與軌道脫離的狀態(tài);輪軌沖擊產(chǎn)生的最大輪軌力為527.5 kN，比軌道平順時增大了2.36 倍，鋼軌和車輪的豎向加速度峰值比軌道平順時分別增大500 倍和40 000 倍，說明軌道不平順是車輛－軌道－橋梁耦合振動系統(tǒng)的重要激勵源之一。

4.2 車輛速度對系統(tǒng)動力響應(yīng)的影響

在軌道不平順相同的前提下，車輛速度在20 m/s ～80 m/s 之間取值，不同車速條件下耦合系統(tǒng)的動力響應(yīng)結(jié)果如表3所示。表3 顯示，車輛速度對輪軌作用力、鋼軌加速度和支承塊加速度的影響很大，且在高速時影響突出;但對車體和橋梁加速度及系統(tǒng)位移幾乎沒有影響，這是由耦合系統(tǒng)傳遞函數(shù)本身性質(zhì)決定的。對于車輛－軌道－橋梁非線性耦合系統(tǒng)來說，車輛速度作為輸入信號主要通過改變輪軌沖擊力，進而改變鋼軌和車輪的加速度對系統(tǒng)起作用，因此，速度對輪軌力響應(yīng)和鋼軌加速度響應(yīng)的影響很明顯。

表3 系統(tǒng)動力響應(yīng)隨車速的變化值

經(jīng)過一系、二系懸掛和軌道二級減振系統(tǒng)的削弱作用，傳遞到車體和橋梁的加速度很小，變化不明顯。直接積分法中位移變化量與加速度變化量之間為時間步長的平方倍關(guān)系，在時間步長足夠小時，輪軌高頻沖擊雖然產(chǎn)生很大的輪軌力和鋼軌瞬態(tài)加速度，但位移變化量幾乎為0。因此，系統(tǒng)位移響應(yīng)受車輛速度的影響很小。

4.3 扣件和膠墊豎直剛度對系統(tǒng)動力響應(yīng)的影響

在軌道不平順相同的條件下，車速取50 m/s，扣件和支承塊下膠墊豎直剛度在20 ～100 kN/mm之間取值，不同剛度條件下耦合系統(tǒng)的動力響應(yīng)如圖9所示。

圖9 鋼軌中點、中間支承塊和橋梁跨中點豎直位移和加速度隨剛度的變化

由圖9 可知，當(dāng)扣件豎直剛度增大時，鋼軌豎直位移減小，最大可減小52%，支承塊位移、軌道和橋梁的豎直加速度都增大而橋梁位移基本不變。這說明扣件剛度不宜過大，否則支撐塊的運動受到橋梁的阻礙，與橋梁發(fā)生激烈的碰撞而導(dǎo)致支撐塊受壓破壞;當(dāng)支承塊下膠墊豎直剛度增大時，軌道和橋梁的豎直位移都減小，分別可減小53%、75%和73%，軌道豎直加速度先增大后減小，橋梁加速度一直增大。這說明膠墊剛度也不宜過大，否則會增大橋梁振動，造成橋梁的損壞。

從整體趨勢來看，系統(tǒng)部件的位移隨剛度增大呈減小的趨勢，振動加速度隨剛度增大呈增大的趨勢。因此，為了取得位移控制和振動控制的平衡點，可以取中間剛度取值，即在40 ～80 kN/mm 之間取值。

5 結(jié)論

本文建立了車輛－軌道－橋梁耦合系統(tǒng)的豎直振動矩陣方程，編寫了計算車輛－軌道－橋梁豎直耦合振動簡化模型的MATLAB 程序，并對程序的有效性進行了算例驗證，其結(jié)論如下:

1)橋梁的振動頻率分布在0 ～200 Hz 的低頻段，以一階彎曲振動為主;軌道的振動頻率分布在0～500 Hz 的中低頻段，以低于20 Hz 的振動為主。

2)在軌道不平順的激勵下，輪軌沖擊力可達(dá)軸重的3 倍，是造成軌道結(jié)構(gòu)振動和鋼軌損傷破壞的重要原因之一;軌道架設(shè)時要對基礎(chǔ)進行夯實處理，防止地基沉降變形引起的軌道幾何形位的變化，以保證車輛安全、平穩(wěn)運行。

3)車輛速度的變化對輪軌沖擊力和軌道加速度的影響很大，對結(jié)構(gòu)位移的影響很小，高速條件下軌道結(jié)構(gòu)的振動控制和安全維護是重要的研究方向。

4)橋上鋪設(shè)彈性支承塊式無砟軌道時，為減小軌道和橋梁的振動，防止支承塊破壞，扣件和支承塊下膠墊豎直剛度要相互匹配，剛度宜在40 ～80 kN/mm 之間取值。

［1］張愛蓮，陳書劍. ADAMS 柔性體建模技術(shù)研究［J］. 煤礦機械，2011(6):95.

［2］晉智斌，強士中，李小珍.高速列車－橋梁豎向隨機振動的時域分析方法［J］.地震工程與工程振動，2008(3):110.

［3］張志超，張亞輝，林家浩.車橋耦合系統(tǒng)非平穩(wěn)隨機振動分析的虛擬激勵－精細(xì)積分法［J］.工程力學(xué)，2008(11):197.

［4］逄煥平，董滿生，侯超群. 車橋耦合振動分析的狀態(tài)空間法［J］.合肥工業(yè)大學(xué)學(xué)報:自然科學(xué)版，2012(12):1610.

［5］羅文俊，雷曉燕，練松良.車輛－高架橋耦合系統(tǒng)豎向振動分析車輛軌道新模型［J］.華東交通大學(xué)學(xué)報，2013(2):1.

［6］李東平，王寧波，曾慶元.車橋時變系統(tǒng)耦合振動分析模態(tài)綜合法［J］.鐵道科學(xué)與工程學(xué)報，2009(3):36.

［7］Vo T P，Lee Jaehong，Ahn N. On sixfold coupled vibrations of thin－walled composite box beams［J］. Composite Structures，2009(89):524.

［8］Nallasivam K，Dutta A，Talukdar S. Dynamic analysis of horizontally curved thin－walled box－girder bridge due to moving vehicle［J］.Shock and Vibration，2007(14):229.

［9］Wu S Q，Law S S. Dynamic analysis of bridge-vehicle system with uncertainties based on the finite element model［J］. Probabilistic Engineering Mechanics，2010(25):425.

［10］翟婉明.車輛軌道耦合動力學(xué)［M］. 北京:中國鐵道出版社，1997.

［11］王寧波，任偉新，肖祥.列車－橋梁耦合振動研究綜述［J］.力學(xué)進展，2012(5):634.

［12］Uzzal R U A，Ahmed A K W，Rakheja S. Analysis of pitch plane railway vehicle-track interactions due to single and multiple wheel flats［J］.Proceedings of the Institution of Mechanical Engineers，2009(1):375.

［13］杜憲亭，夏禾，張?zhí)?車橋耦合振動迭代求解穩(wěn)定性研究［J］.振動與沖擊，2012(22):62.

［14］Zhua J J，Ahmed K W，Rakhej S，et al. Development of a vehicletrack model assembly and numerical method for simulation of wheel-rail dynamic interaction due to unsupported sleepers［J］.Vehicle System Dynamics: International Journal of Vehicle Mechanics and Mobility，2010(12):1535.

［15］李永樂.風(fēng)－車－橋系統(tǒng)非線性空間耦合振動研究［D］. 成都:西南交通大學(xué)，2003.

［16］Lou Ping，Yu Zhiwu，Au F T K. Rail-bridge coupling element of unequal lengths for analyzing train-track-bridge interaction systems［J］.Applied Mathematical Modelling，2012(36):1395.

［17］李小珍，強士中.列車－橋梁耦合振動研究的現(xiàn)狀與發(fā)展趨勢［J］.鐵道學(xué)報，2002(5):112.

［18］徐榮橋.結(jié)構(gòu)分析的有限元法與 MATLAB 程序設(shè)計［M］. 北京:人民交通出版社，2005.