基于PFC-FLAC耦合的彈性軌枕力學特性分析

崔旭浩，肖 宏

(北京交通大學 軌道工程北京市重點實驗室, 北京 100044)

近年來，隨著列車運行速度和軸重的不斷提高，輪軌相互作用不斷加劇，散體道床在劇烈振動下越發容易出現道砟破碎、粉化等劣化現象，這不僅增大了工務養護維修工作量，還會產生道床板結等影響行車安全，因此采取有效的措施減小道砟破碎、粉化是十分必要的。彈性軌枕作為一種保護道床的有效措施，在國內外均有一定應用。針對彈性軌枕開展力學特性研究非常有必要。

在試驗研究方面，馬春生等[1]進行了高速鐵路橋上彈性軌枕有砟道床靜、動力學特性試驗研究，其認為在高速鐵路橋上采用彈性軌枕道床可以滿足高速鐵路的開通要求。練松良等[2]開展了軌道結構的動力性能試驗，其認為枕下膠墊可以有效的調整軌道結構的整體剛度。Omodaka等[3]針對鋪設彈性軌枕后道床沉降、橫向阻力以及枕下膠墊的耐久性進行了試驗研究。Sol-Snchez等[4]以廢舊輪胎制成的減振墊為對象，通過試驗探討了其在軌枕下鋪設的可行性，其研究結果表明通過改變枕下膠墊的厚度及材料的加工工藝可以獲得不同靜態和動態剛度的彈性墊。Kaewunruen等[5]研究了鋼軌接頭處鋪設彈性軌枕對軌道結構的影響，結果表明在鋼軌接頭處使用彈性軌枕會加大鋼軌和軌枕的振動。可以看出，針對彈性軌枕的試驗研究主要集中在軌道結構的靜、動態響應及材料的適用性等方面，針對彈性軌枕與碎石道砟的相互作用方面幾乎沒有涉及。

在理論研究方面，Johansson等[6]通過數值模擬研究分析了枕下膠墊對輪軌動力相互作用的影響，其研究表明影響主要在0～250 Hz的范圍內。蔡成標等[7]通過模擬落軸試驗，探討了軌下剛度和塊下剛度的合理取值問題。焦雷[8]基于車輛-軌道-下部基礎空間耦合動力學，運用有限元軟件Abaqus從行車安全性、平順性和軌道減振角度，對彈性長枕減振軌道結構的減振特性進行了研究。翟淼[9]根據車輛-軌道耦合動力學理論，運用Ansys對橋上有砟軌道鋪設彈性軌枕的豎向振動特性進行了研究。Cao等[10]建立了車輛-軌道的半空間模型，對臨界行車速度下的彈性軌枕對地面的隔振作用進行了分析。Gao等[11]基于離散單元法，從道床沉降、道砟振動等角度出發，運用PFC3D研究分析了鋪設彈性軌枕后有砟道床服役性能的改善狀況，但其模型中未考慮軌枕的變形特性。可以看出，既有針對有砟道床鋪設彈性軌枕的理論研究方面，主要可以劃分為兩種：①基于有限單元法，將散體道床簡化為連續的彈性體，這種方法可以較好的模擬軌枕，但對于散粒體結構的有砟道床模擬則存在缺陷，難以考慮道砟顆粒復雜外形及道床的級配特征；②基于離散單元法，這種方法可以較好的模擬有砟道床的散體特性，但對于軌枕這一連續體則只能簡單的視為邊界墻單元或塊體進行模擬，其不能反映軌枕的受力變形特征。此外，從上述的理論研究現狀還可看出，對彈性軌枕有砟道床的研究多是從軌枕及道床的振動指標及沉降特性角度進行，針對彈性軌枕及道砟顆粒的受力特征，例如彈性墊層與道砟顆粒之間的接觸應力、彈性軌枕對有砟道床內部道砟顆粒之間接觸力的影響等方面的研究基本空白，而軌枕和道床的受力直接影響其服役狀態。

基于此，本文將針對既有研究的不足，率先采用離散元(Particle Flow Code，PFC)與有限差分法(FLAC)聯合仿真的方法建立彈性軌枕-散體道床相互作用模型，對有砟道床鋪設彈性軌枕的力學特性開展深入研究。

1 計算模型的建立

1.1 有砟道床離散元模型的建立

采用PFC進行道砟顆粒的建模。PFC是離散單元法的一種，用于顆粒材料在二維和三維條件下的力學分析，其基本單元是圓盤和圓球顆粒，利用邊界墻(Wall)進行約束。PFC中計算時不需要給材料參數定義宏觀本構關系和對應的參數，而是采用局部接觸反映宏觀問題，可以研究分析顆粒間的相互作用問題、大變形問題、斷裂問題等。

精確模擬道砟顆粒的外形是有效分析有砟道床力學特性的關鍵所在。趙春發等[12]通過對碎石道砟樣本照片進行圖像處理，并采用六邊形密排圓盤填充的方法，實現了道砟的模擬，但得到的道砟數值模型內部密排的圓盤之間存在孔隙，這與現實中致密的道砟材質不相符。Xiao等[13]采用三維激光掃描的方法對道砟顆粒的三維外形信息進行了獲取，通過投影的方法得到了道砟顆粒的二維輪廓，并建立了道砟顆粒的離散元數值模型，但在道砟顆粒模擬時人為引入孔隙的問題依舊存在。針對這一問題，本文使用數碼相機對碎石道砟樣本進行拍照，運用Matlab進行圖像處理獲得道砟顆粒的外輪廓信息并導入到PFC2D中，以導入的道砟顆粒外輪廓為邊界通過基本圓盤單元(Pebble)堆疊的方式生成組合道砟顆粒的塊體模型。建模過程見圖1。

圖1 道砟顆粒建模過程

通過參數設置改變組合道砟的基本單元個數可以實現不同精度道砟外形的模擬，不同精度道砟顆粒的模擬情況見圖2。由圖2可知，組成塊體的基本單元數量越多則道砟計算模型越逼近真實的道砟外形，但計算耗時也會大幅增加。經過多次試算，權衡模擬精度和計算耗時后，選用32個基本單元進行堆疊的方式模擬組合道砟顆粒。在仿真計算中依據生成道砟塊體的輪廓邊界和道砟密度計算道砟顆粒的質量、轉動慣量及質心位置。根據TB/T 2140—2008《鐵路碎石道砟》[14]，采用特級道砟級配建立模型，分別見圖3、圖4。參考文獻[15-17]，取道床模型參數見表1。

圖2 不同精度道砟計算模型

圖3 特級道砟粒徑級配曲線

圖4 有砟道床模型

表1 離散元模型參數

1.2 軌枕有限差分模型的建立

彈性軌枕由“普通軌枕+枕下墊層”組成，建立的有限差分FLAC2D彈性軌枕模型，見圖5，有限差分模型計算參數見表2。為更好的分析彈性軌枕的力學特性，同時建立普通軌枕-道床模型進行對比分析。

圖5 彈性軌枕模型

表2 有限差分模型參數表

1.3 離散元與有限差分的耦合

周健等[18]通過PFC和FLAC的Socket I/O接口，針對基本圓盤單元與連續介質單元的接觸位移關系進行了分析，最終實現二維離散-連續的耦合計算，但其并未考慮顆粒具有不規則外形的情況。本文將在論述組合道砟顆粒與連續單元交界面接觸位移計算原理的基礎上，通過Socket I/O接口實現PFC和FLAC的耦合計算。

組合道砟顆粒與連續單元的接觸示意見圖6。

圖6 組合道砟顆粒與連續單元接觸示意

組合道砟顆粒的基本單元與交界面單元的接觸點坐標為

(1)

式中：xBi為組合道砟顆粒基本單元的形心坐標；di為基本單元中心到交界面單元平面的垂直距離；RBi為基本單元的半徑；ni為交界面單元的內法線法向。

組合道砟顆粒與交界面單元的接觸作用力FC為

FC=Fn+Fs

(2)

式中：Fn為法向分量；Fs為切向分量。

法向分量計算公式為

(3)

式中：m為組合道砟顆粒中與交界面單元相接觸的基本單元個數；Ui為基本單元與交界面單元的“重疊”量，Ui=RBi-di；Kn為接觸點的法向剛度。

切向分量以增量形式計算，組合道砟顆粒基本單元與交界面單元在接觸點的相對運動速度為

(4)

在每一時步Δt內接觸位移增量矢量為

Δxi=viΔt

(5)

接觸點的切向位移增量分量為

Δxsi=Δxi-Δxini

(6)

組合道砟顆粒在交界面上所受切向接觸力增量矢量為

(7)

式中：ΔFs為切向接觸力增量；Ks為接觸點的切向剛度。接觸形成時總的切向接觸力為零，通過累計切向接觸力增量計算切向接觸力為

(8)

組合道砟顆粒所受接觸合力矩M計算公式為

(9)

在獲得組合顆粒與交界面單元的接觸力后，通過等效節點荷載分配的方法計算交界面單元的節點荷載。分別到離散單元體系和連續單元體系中根據運動平衡方程計算速度位移，在下一步分析開始時，根據上一步的計算結果更新交界面單元和組合道砟顆粒的位置，重新計算接觸作用，再回到各自系統中進行這一時步的運算。如此迭代，實現連續時間的仿真計算，計算流程框架見圖7。參考文獻[9]，確定耦合邊界的參數取值見表3。

圖7 耦合邊界計算流程框架

表3 耦合邊界參數

2 模型的驗證及外荷載的施加

2.1 模型的驗證

為了驗證模型的正確性，選定鋪設彈性軌枕及普通軌枕的有砟道床試驗段開展動態測試，測試過程中列車運行速度為60 km/h。實測鋪設彈性軌枕的道床振動加速度峰值為5.68g，實測普通軌枕道床的振動加速度峰值為6.63g。將現場實測的枕上動壓力時程作為輸入荷載施加到本文所建模型中，計算得彈性軌枕道床振動加速度峰值為5.57g，普通軌枕道床振動加速度峰值為6.56g。現場實測的彈性軌枕道床振動加速度時程圖和仿真計算的彈性軌枕道床振動加速度時程圖，見圖8。由圖8可知，實測和仿真計算結果波形相似，加速度峰值也十分接近，說明了本文所建立模型的正確性。

圖8 道床加速度實測結果和計算結果

2.2 外荷載的施加

本模型計算時在軌枕承軌槽位置處施加垂向簡諧荷載為

F(t)=F0[sin(2πft-π/2)+1]

(10)

式中：F0為荷載半峰值；f為荷載頻率；t為荷載作用時間。參考文獻[19]，本文選取荷載半峰值為280 kN，將荷載均分施加在軌枕頂面的承軌槽處。對于2.4 m固定軸距貨車，行車速度為60 km/h時，荷載作用頻率的理論值約為7 Hz，故本文施加的動力荷載的頻率取為7 Hz。

3 計算結果分析

3.1 軌枕的受力特性

荷載加載至峰值時普通軌枕和彈性軌枕的垂向動應力分布云圖，見圖9。

圖9 不同軌枕垂向動應力分布

由圖9可知，兩種軌枕底部因其與散體道砟的非連續接觸而導致應力都呈現非均勻分布。在相同的接觸位置處，彈性軌枕底部的動應力幅值明顯小于普通軌枕，這主要是由于彈性軌枕的枕下墊剛度小，在與道砟相接觸的位置在受力后會發生變形、增加接觸面積，從而減小接觸應力。彈性軌枕與散體道砟的細觀接觸見圖10，由圖10可知，在與道砟相接觸的位置，彈性墊層發生了明顯的變形。

圖10 耦合邊界細觀圖

提取軌枕底面2.6 m范圍內與散體道砟的接觸應力見圖11。由圖11可知，普通軌枕與彈性軌枕底部接觸應力的分布呈現出相同的趨勢，在線路中心線位置處的應力較小，鋼軌對應位置處應力較大，彈性軌枕與道床的接觸應力幅值為228 kPa，普通軌枕底面的接觸應力幅值為435 kPa，彈性軌枕底面的接觸應力幅值比普通軌枕小47.5%。軌枕對道床的接觸應力越小，道砟越不易破碎，道床的養護維修工作量就越小。由此可見，相對于普通軌枕，彈性軌枕對道床有良好的保護作用，能大大減小工務養護維修工作量。

圖11 軌枕底面與道砟顆粒接觸應力分布

為分析軌枕在荷載作用下的受力特性，提取荷載作用下軌枕頂面、中上部(距離上表面0.06 m)、中部(距離上表面0.12 m)、中下部(距離上表面0.18 m)、底面位置處的垂向動應力的分布特征。荷載施加至峰值時軌枕內部各位置處的垂向動應力沿橫向的分布情況，見圖12。

圖12 軌枕垂向動應力沿橫向分布圖

由圖12可知，峰值荷載作用下，彈性軌枕和普通軌枕內部不同位置處動應力沿橫向分布的規律基本一致，軌枕頂部承軌槽位置處所受垂向動應力最大。該規律與既有認識一致，從而表明采用FLAC與PFC進行有砟軌道結構分析的方法正確性。

3.2 道床的受力特性

有砟道床內部道砟顆粒之間的接觸力大小會影響其劣化及壽命，過大的接觸力會加快道砟顆粒的磨耗及破碎的發生。

循環荷載加載至峰值時普通軌枕道床和彈性軌枕道床道砟顆粒接觸力的分布情況，見圖13。

圖13 不同有砟道床力鍵分布

由圖13可知，兩種道床結構中強力鍵都主要分布在軌枕的正下方。普通軌枕道床內強力鍵在道床頂層范圍內分布比較集中，尤其是在軌枕相接觸的位置，而彈性軌枕道床頂層范圍內強力鍵出現的個數明顯較少，這表明彈性軌枕可以減小道床頂層范圍內道砟顆粒之間強力鍵的個數。在道床底層范圍內，兩種道床結構中強力鍵的分布情況差異較小。

統計加載之前和荷載峰值作用時道砟顆粒的平均和最大法向、切向接觸力，見表4。

表4 道床接觸力統計表 kN

由表4可知，僅在加載前自重作用下，道砟顆粒之間的接觸力普遍較小，彈性軌枕道床中的各個接觸力指標與普通軌枕道床基本相等，而道砟顆粒法向接觸力明顯大于切向接觸力，前者約為后者的3倍，說明自重作用下，道床內部的荷載傳遞主要以法向接觸傳力為主。

當荷載加載至峰值時，道砟顆粒之間的接觸力較自重作用下有大幅度提高，約提高了一個量級。對比數據可知，彈性軌枕道床道砟顆粒平均法向接觸力和最大法向接觸力分別比普通軌枕道床減小7.0%和13.4%。彈性軌枕道床道砟顆粒平均切向接觸力和最大切向接觸力分別比普通軌枕道床減小8.8%和9.3%。由此表明，鋪設彈性軌枕可以顯著降低道砟顆粒之間的平均接觸力和最大接觸力。這主要是由于彈性軌枕枕下彈性墊在動荷載作用時會發生變形，增大道砟與軌枕的接觸面積，從而可以減緩軌枕與道床之間的相互作用力。

在細觀角度評價顆粒之間接觸力方面，有學者提出采用概率密度的方式來度量顆粒間力鍵的大小及分布情況[20-21]。接觸力的概率密度可以衡量顆粒之間的接觸力出現在某一范圍的概率。本文對道砟顆粒之間的接觸對進行遍歷，統計接觸對的個數和接觸力的大小，力鍵的概率密度為

(11)

(12)

式中：x為道砟顆粒之間的接觸力，實際計算中x取不到+∞，而是以道砟顆粒之間的最大接觸力為取值區間上限；ν(x)為所有接觸對中接觸力小于x的接觸對的個數；n為接觸對總數；F(x)為概率分布函數；p(x)為概率密度函數。

自重作用和荷載加載至峰值時道床中法向接觸力和切向接觸力的概率密度分布，見圖14。由圖14可知，不論是法向接觸力還是切向接觸力，隨著顆粒間接觸力增加，其概率密度逐漸減小，道床內部出現較大接觸力的概率較小，這表明在道床內部大部分道砟顆粒之間的接觸力較小。

圖14 道砟接觸力的概率密度分布

由圖14可知，當荷載加載至峰值時，在彈性軌枕道床和普通軌枕道床中，弱力鍵的概率密度數值十分接近，但強力鍵在兩種道床結構中的概率密度分布有明顯差異，表現為當顆粒間法向接觸力大于7.626 kN時，彈性軌枕道床中的法向接觸力概率密度開始低于普通軌枕。切向接觸力亦有同樣類似的規律，峰值荷載作用下，當接觸力超過2.641 kN時，彈性軌枕道砟顆粒切向接觸力的概率密度開始小于普通軌枕。由此表明，鋪設彈性軌枕可以減低道床內部較大接觸力出現的概率，但該現象只有在荷載達到一定值時才能體現出來。由此可見，若將彈性軌枕鋪設在重載鐵路上，對道砟的保護作用將更加顯著。

3.3 振動加速度響應分析

為探究彈性軌枕對道床結構振動響應的影響，對動荷載作用下軌枕和道床的振動加速度進行計算分析，測點布置見圖15。

圖15 振動加速度測點布置(單位：mm)

繪制出計算位置處軌枕和道床振動加速度的時程曲線，分別見圖16和圖17。

圖16 軌枕振動加速度

圖17 道床振動加速度

由圖16和圖17可知，兩種道床結構中軌枕和道床振動加速度時程圖波形相似，但兩者的振動加速度峰值存在差異，普通軌枕的振動加速度峰值為1.21g，彈性軌枕的振動加速度峰值為1.43g；普通軌枕道床的振動加速度峰值為4.95g，彈性軌枕道床的振動加速度峰值為2.76g。對比可知，彈性軌枕振動加速度比普通軌枕大18.2%，彈性軌枕道床振動加速度比普通軌枕道床小44.2%。由此可見，鋪設彈性軌枕會增大軌枕的振動加速度而減小道床的振動加速度，且對道床的振動加速度影響更為顯著。另外，還可以看出，道床的振動加速度大于軌枕的振動加速度，該規律與文獻[12]的研究規律一致，這是由于道砟的參振質量較小，更易引起強烈振動。

3.4 彈性墊板合理靜剛度的探討

在彈性軌枕設計時，枕下膠墊的靜剛度是一個重要參數，其直接關系到彈性軌枕-有砟道床的工作狀態。為分析枕下膠墊靜剛度對道床動態響應的影響，通過改變模型中膠墊耦合邊界的接觸剛度參數，分別取靜剛度為10、20、30、…、90、100、…、290、300 kN/mm的彈性膠墊進行參數影響分析。參考文獻[9]，按照線性關系對軌枕模型中彈性膠墊的楊氏模量參數進行調整。軌枕振動加速度、道床振動加速度和道床彈性變形隨枕下膠墊靜剛度的變化情況見圖18～圖20。

圖18 軌枕振動加速度

圖19 道床振動加速度

圖20 道床彈性變形

由圖18可知，隨著枕下膠墊靜剛度的逐漸增大，軌枕的振動加速度逐漸減小，但當膠墊的靜剛度大于40k N/mm時，軌枕振動加速度下降變緩，并保持一個較低的范圍、趨于收斂。可見，若要使軌枕的振動在較低水平，膠墊的靜剛度應大于40 kN/mm。

由圖19可知，當膠墊的靜剛度小于100 kN/mm時，道床的振動加速度穩定在較低范圍內；當膠墊的靜剛度超過100 kN/mm時，隨著膠墊靜剛度的增大，道床的振動加速度逐漸增大。因此，為使道床的振動不過于劇烈，建議膠墊的靜剛度小于100 kN/mm。

由圖20可知，隨著枕下膠墊靜剛度的逐漸增大，道床的彈性變形逐漸減小。當膠墊的靜剛度為70 kN/mm時，達到了TB 10761—2013《高速鐵路工程動態驗收技術規范》[22]中規定的基準值1.5 mm。因此，從控制道床彈性變形的角度，應使得膠墊的靜剛度大于70 kN/mm。

綜上，建議枕下膠墊的靜剛度取70～100 kN/mm。

4 結論

(1)本文采用離散單元法建立了具有棱角特征的精細化道床模型，采用有限差分法建立了可考慮變形的連續介質軌枕模型，構建了組合道砟顆粒與有限差分單元之間的接觸位移關系，通過在軌枕與道床之間設置耦合邊界進行力和位移的交互實現二者的耦合計算，該模型能很好的進行軌枕與有砟道床的相互作用分析。

(2)軌枕底面與散體道砟之間存在非連續接觸從而導致軌枕底面應力呈現非均勻分布，軌枕底面在與道砟相接觸的位置會出現較大的接觸應力。但由于彈性軌枕枕下彈性墊的剛度較小，受力后會發生變形、增加接觸面積，從而可以減小軌枕與道砟之間的接觸應力。

(3)彈性軌枕可以減少道床頂層范圍內道砟顆粒之間較大接觸力，對底層道砟顆粒的接觸力影響較小。彈性軌枕可以減小道砟顆粒之間的平均接觸力和最大接觸力，可以起到保護道砟顆粒的作用。

(4)在道床內部，道砟顆粒之間的接觸力的概率密度隨著接觸力的增加而減小，在道床內部大部分道砟顆粒之間的接觸力較小。另外，鋪設彈性軌枕可以減低道床內部較大接觸力出現的概率，在荷載較大時更為明顯，從減小接觸力避免道砟破碎角度，建議重載鐵路在必要時使用彈性軌枕。

(5)彈性軌枕可以顯著降低道床的振動加速度，同時軌枕本身的加速度稍有增大。從合理控制軌枕、道床的振動加速度和彈性變形綜合考慮，建議枕下膠墊的合理靜剛度宜為70～100 kN/mm。