基于離散元與有限元耦合的有砟道床沉降特性研究

王 桐，王立華，陳佳明，栗先增

(昆明理工大學 機電工程學院，云南 昆明 650500)

鐵路作為我國的重要基礎設施，在交通運輸體系中占據重要地位，也是目前應用最廣泛的運載方式之一，對我國的經濟建設和發展有著十分重大的意義。由于我國鐵路近幾年發展迅速，且與各個領域之間的交叉融合越來越多，所以對于鐵路方面的研究也愈加深入。有砟鐵路作為我國鐵路運載的重要組成部分，是鐵路發展的基礎，對其深入的研究就顯得尤為重要。有砟軌道道床質量對于保障列車運行穩定安全方面具有重要意義。

在對有砟鐵路進行研究時，由于組成道床的道砟屬于離散體顆粒，而道床上部軌枕屬于連續體，所以在對二者進行分析時不能使用單一的方法。有限元法只能從宏觀上反應軌枕與道床的整體特性。Ishikawa[1]等通過室內試驗研究有砟道床在定點載荷和移動載荷作用下的力學行為。其結果表明，在移動載荷作用下有砟道床產生的塑性變形比定點位移作用下的大。Mc Dowell[2-4]等基于真實顆粒外形的最小半徑離散單元數值模型的構造方法，用所生成的組合球單元對復雜外形的道砟顆粒進行了模擬，并對比分析了組合單元顆粒與球形單元顆粒在荷載作用下的力學特性。井國慶[5]等基于二維離散單元法顆粒軟件Particle Flow Code 2D(PFC 2D)，建立了道砟軌枕離散單元模型，研究了列車循環荷載作用下道砟顆粒破碎與軌枕沉降的關系。Shaer[6]等采用等比例縮減尺寸的鐵軌-枕木-道床-路基模型得到整體道床的沉降與枕木和道砟的振動加速度密切相關。高亮等[7-8]利用離散單元法，分別建立球形道砟顆粒和顆粒簇形成道砟顆粒的軌枕與散體道床的離散元模型，對比分析了兩種不同模型在循環荷載下的力學性能，并利用離散元顆粒數值模型構造道砟顆粒模型。同時，在此基礎上建立了循環荷載道砟箱數值模型，研究在列車高速及重載線路條件下，道砟級配對散體道床沉降力學特性的影響情況，并從細觀角度分析了散體道床的沉降原理。但若要同時從宏觀和微觀方面對有砟道床進行分析，就需要引入有限元和離散元耦合的方法。

針對離散元與有限元耦合方法的研究，馮春[9]將對有限元模型進行加載達到臨界狀態時，把有限元替換為離散元研究。嚴穎等[10-12]在離散元與有限元的接觸面處提出一種新的接觸算法，實現有砟鐵路軌道的離散元-有限元耦合法的研究分析，又通過鑲嵌單元模擬來達到耦合的方式并通過Fortran語言編程計算有限元與離散元的耦合結果。嚴穎等依據接觸面上的能量守恒原理，確定兩種介質間的接觸力實現離散元-有限元單元區間力子參數傳遞。尹超[13]利用MATLAB編制程序接口，通過間接計算實現離散元軟件與有限元軟件的數據計算與傳遞。

上述文獻都是基于耦合方法研究軌枕與散體道砟顆粒的相互作用力，需要通過外部程序來進行有限元和離散元顆粒之間的數據傳遞。為簡化數據之間傳遞的繁瑣過程，本文采用有限元軟件得到軌枕網格坐標，導入離散元中用離散元顆粒以及鍵連接構造軌枕模型。后采用球顆粒組合模型構造非規則道砟，得到縮小比例工程尺度的枕木-道床離散元模型。本研究通過對軌枕施加循環正弦激振，研究了軌枕和散體道砟顆粒的相互作用和軌枕的受力、變形情況。此外，本文還研究了不同激振頻率對有砟道床力學性質的影響，得到道床的沉降規律并探究了道床累積沉降的內在機理。

1 有砟道床的離散元模型

1.1 典型形狀道砟模型構建

離散單元法(Discrete Element Method，DEM)主要用于離散顆粒物料研究領域。非規則形狀顆粒建模耗費大量時間，且模型精度對顆粒間的動力學響應特性有較大影響，故本文采用離散單元法構建道床模型，旨在提高非規則形狀顆粒建模效率及精度。在有限元軟件中建立軌枕模型并劃分網格，同時獲取網格單元坐標信息，并在離散元軟件中設置接觸半徑及物理半徑，然后通過工程離散單元法(Engineering Discrete Element Method，EDEM)中的應用程序編程接口實現多球快速聚合顆粒建模。道砟顆粒的形態會對散體道床的動力特性產生顯著的影響，構成簇顆粒的球單元越多，其對道砟形態的描述就越真實。利用以上構造方法可以對任意形態的道砟進行構造。考慮到離散元計算效率，本文采用4～10顆小球顆粒組成的簇顆粒單元對不規則道砟顆粒進行近似構造，模型示意圖如圖1所示。組成簇顆粒的球顆粒之間存在很大的重疊量，且保證了道砟內部不存在孔隙。顆粒簇在離散元模擬計算中不考慮其破碎，每個顆粒簇在運動過程中被視為一個整體，顆粒簇之間采用線性接觸模型去模擬顆粒之間的接觸行為。

圖1 不同形狀的鐵路道砟顆粒簇模型

正常工作的有砟鐵路使用的道砟顆粒尺寸一般在37～54 mm之間，如上圖所示粒徑范圍分別在表1道砟級配中的30～45 mm、45～60 mm、>60 mm這3個范圍內選取。

表1 道砟級配

道砟與道砟之間，道砟與構成軌枕的離散元顆粒間的接觸參數如表2所示

表2 道砟間的接觸參數

1.2 節點坐標提取與離散元軟件顆粒的替換

在離散元中將分析對象視為充分數量的離散介質進行計算，其將每個離散顆粒作為一個單元，然后求解出每一時刻中顆粒的相互作用、接觸力及離散單元的運動狀態。隨后，通過引入時間步長的方式，依次求解出離散單元體的運動形態。離散單元體一般使用球形顆粒，根據實際問題的不同，將顆粒模型分為軟球模型和硬球模型兩種，本文所采用的是軟球模型。

利用離散元軟件EDEM所提供的應用程序編譯外接口(Application Programing Interface，API)編寫插件來自定義接觸模型、顆粒生成、外部耦合及顆粒模型[14]。軌枕顆粒建模仿真實驗中利用Visual Studio編程軟件創建的顆粒工廠動態鏈接庫文件，得到顆粒ID、顆粒三維坐標、顆粒實際物理半徑，顆粒接觸半徑等所需的模板顆粒信息，并對顆粒進行替換。考慮到顆粒替換得到的模型與顆粒之間的粘結鍵及仿真所需要的計算時間，文中將組成軌枕的顆粒小球實際物理半徑設置為1.25 cm，顆粒接觸半徑設置為1.6 cm。由離散元顆粒替換構成的軌枕模型如圖2所示。

圖2 離散元顆粒替換構成的軌枕

1.3 軌枕-散體道床耦合模型的建立

將道砟顆粒模型建模后，根據已有相關研究[15]和道砟顆粒模型的尺寸建立道砟箱模型。當道砟容器的尺寸大于顆粒平均粒徑的8倍時，容器的邊界效應可以被忽略。本文建立的道砟顆粒模型的平均粒徑不超過45 mm，容器最小要求尺寸為340 mm，建立鐵路道砟箱三維模型，其中道砟箱的尺寸為1 000 mm×700 mm×600 mm(長×寬×高)。考慮到仿真所面臨的時間和內存問題，本文僅建立由一根軌枕和有砟道床組成的離散元模型。簡化模型中有砟道床及軌枕的尺寸與實際道床及軌枕的尺寸比例為1∶8，軌枕由簡單的長方體簡化代替。

在有限元軟件中對節點數據進行處理后，在離散元軟件中，設置單獨母體小球，在到達需要放置的軌枕位置后，通過顆粒替換將離散元顆粒按照軌枕網格劃分坐標在離散元軟件中排布組合成新軌枕。建立符合要求的有砟道床離散元模型，如圖3所示。

圖3 軌枕與道床的三維離散元模型

完成軌枕與道床的離散元模型建立后，為其添加接觸參數，確定軌枕顆粒和道砟顆粒運動邊界條件。由于對道砟顆粒的研究需要確定顆粒間準確的接觸力變化情況，且暫不考慮顆粒破碎的情況，因此采用Hertz-Mindlin(no slip)無滑動接觸模型。由于軌枕由散體顆粒構成，需要設置bond鍵來連接各個離散元顆粒以達到整體作為軌枕的效果。軌枕采用Hertz-Mindlin with bonding基礎模型。用于模擬破碎、斷裂等問題時，則采用小顆粒粘結成大塊物料。外力作用下顆粒間粘結力會發生破壞，從而產生破碎及斷裂效果，但本文不考慮軌枕的劣化和破損問題，故bond鍵強度設置較大。最終所確定的道砟顆粒與顆粒，顆粒與軌枕間的恢復系數、靜摩擦系數和動摩擦系數如表3所示。

表3 軌枕粘結模型物理參數

由于無法直接對軌枕施加激振載荷，因此需要將作用在軌枕上的載荷轉化為不同的形式施加在軌枕上。將軌枕所受到的激振力轉化成位移的形式輸入到軌枕上，即在仿真過程中通過在軌枕兩側施加激振源，對激振源設置位移函數，水平激振力屬于正弦激振力，存在頻率和幅值。因此將水平激振力轉化成頻率和位移幅值的形式施加在軌枕上，并通過改變激振頻率來對不同激振頻率下軌枕的沉降狀態進行研究。對軌枕施加激振，完成穩定過程仿真，以離散元軟件為主體來實現耦合過程。

2 穩定作用下軌枕受力及道床動力特性分析

2.1 軌枕下方不同區域的速度變化趨勢

在軌枕下方的道砟在穩定作用下，道砟的移動方向整體向下，其中以軌枕下方的道砟顆粒尤為明顯。為了直觀分析不同載荷頻率對道床的沖擊程度，如圖4所示，給出了載荷頻率為36 Hz時枕木受到穩定作用，道床中間厚度1 cm垂向截面區域內所有道砟顆粒的速度分布圖。可看出枕木兩側道砟顆粒的速度均較小，說明穩定載荷對軌枕下方道砟起主要作用。道砟顆粒的速度幅值呈現近似扇形的衰減，且軌枕下方道砟顆粒具有更加明顯的垂直向下運動趨勢，越靠近軌枕中心的道砟顆粒，其垂向速度響應越大。

軌枕在受到激振時，向左右兩側產生微小位移，在軌枕向左向右移動時，砟肩部分的道砟隨著軌枕產生同向位移。由于軌枕下方道砟顆粒始終為向下移動，所以在激振影響中，道砟不斷向下移動，達到重新排布效果。整個過程中，軌枕下方道砟不斷下沉，使上方軌枕下沉。

2.2 軌枕的力學特性分析

選取軌枕底部顆粒，提取該顆粒在穩定過程中的位移量來得到軌枕整體的沉降量，如圖5所示，顆粒在穩定過程中的沉降量為20 mm左右。由圖中可以看出，在對不同激振頻率下道床穩定作業過程的仿真中，隨著穩定作業的進行道床軌枕沉降量逐漸增大，并且隨著激振頻率的增加軌枕總沉降量逐漸增大，在36 Hz時達到最大值；當激振頻率繼續增大時，沉降量又逐漸減小。

圖5 軌枕沉降量與時間的關系

3 設計道砟箱試驗與模型驗證

3.1 道砟箱試驗的設計

根據前文軌枕和道砟顆粒的選型設計搭建試驗平臺，如圖6所示，建立與仿真同比例的道砟箱試驗。文中所使用的箱體尺寸(長×寬×高)為1 000 mm×700 mm×550 mm，所用材料為3 mm厚鋼板。通過安裝偏心振動電機在水泥軌枕上部來模擬穩定過程，并使用調頻器達到不同激振頻率的效果。該振動電機可以提供水平方向的激振力，在垂直下壓力方面只考慮軌枕和振動電機的質量。

圖6 道砟箱試驗

為研究試驗中模擬穩定作業軌枕沉降量的變化規律，通過激光測距儀可以準確地得出動力穩定作業過程中道床的垂向位移和軌枕的水平激振幅值以及水平激振頻率與軌枕沉降的關系。將激光測距儀安裝在軌枕正上方距離3 cm處，同時避免與道砟箱接觸，防止振動對測量產生影響。在進行測量時，測距儀的激光照射到軌枕的振動電機安裝鋼板上，當模擬穩定作業過程開始時，激光測距儀發出紅外線照射到軌枕上方的平滑鋼板上進行數據采集，得到每一個時刻軌枕與激光發射口的距離大小。將初始位置設置為零點，即可得到軌枕的沉降量。

圖7 激光測距儀實物圖

動力穩定車在對劣化道床進行穩定作業的過程中，水平激振頻率的工作頻率范圍是0～45 Hz，劉暢[16]得到的動力穩定作業過程中最優激振頻率為30 Hz；項永志[17]通過仿真得到動力穩定作業最優激振頻率為31～33 Hz。本文建立的試驗平臺為小型道砟箱試驗平臺，由于激振電機的頻率限制，激振頻率在30 Hz以下的激振作用振動比較微弱，所對應的穩定效果不明顯，因此最終選擇設置梯度遞增的激振頻率，分別為30 Hz、32 Hz、34 Hz、36 Hz、38 Hz和40 Hz，共6組。由于偏心振動電機的額定頻率為50 Hz，無法直接改變其頻率，因此通過變頻器來調節電機的頻率，達到本文對于不同頻率的試驗要求。

3.2 數據的提取與驗證

在比例縮小模型實驗中，使用激光測距儀測量軌枕下沉量，與上文中的軌枕下沉量進行對比分析。為方便統計數據和計算，在仿真模擬中選取如圖8所示軌枕底面小球顆粒編號11 316作為分析對象，提取此處顆粒在穩定過程中的位移來進行數據對比分析。

圖8 道砟顆粒選取位置

提取此位置道砟顆粒在穩定過程中的位移即軌枕的沉降量，如圖9所示得到仿真過程中軌枕顆粒的位移情況，可以發現在穩定過程中，軌枕的沉降量隨著激振頻率的增加而增大，當激振頻率為36 Hz時達到最大值，之后隨著頻率的增加逐漸減小。

由圖9可知激振頻率不變時，動力穩定作業下道床沉降趨勢與循環載荷作用產生的道床沉降趨勢基本相同。由于道砟顆粒相互接觸時具有一定的彈性，因此道床沉降過程具有一定的反復，試驗結果和文獻[16]中動力穩定下的道床沉降趨勢基本相同。這是因為道床初期比較松散，密實度低，穩定作業下道砟顆粒之間迅速密實，沉降急劇增大，隨著穩定作業的進行，道床逐漸密實，穩定性逐漸提高，道床沉降量逐漸變緩；當激振頻率增加時，道床沉降量呈現先增大后減小的趨勢，在36 Hz時達到最大，其道床的沉降穩定效果最好。

圖9 道床沉降趨勢圖

將實驗得到的軌枕隨時間變化的沉降量變化情況繪制成軌枕沉降量相對于時間變化關系圖像，并繪制經過仿真得出的軌枕沉降與時間關系圖像對比實驗及仿真結果。如圖10所示，將離散元有限元耦合建模與試驗獲取的數據結果進行對比可知，試驗獲得的軌枕沉降曲線由于在道砟箱堆放道砟時道砟顆粒之間空隙較大堆疊稀松，所以在激振在初期沉降量上升較快。試驗與仿真的沉降趨勢基本一致，說明本文建立的有砟道床仿真分析模型是可靠的，所選取的參數正確，所確定的邊界條件合理，仿真分析結果有效。

圖10 實驗與仿真對比圖

經過對軌枕在36 Hz激振力作用下仿真得到的軌枕沉降量與時間變化可看出其趨勢變化規律形如指數函數。為了方便研究軌枕沉降量變化規律并得到普遍性結論，文中對仿真得到的軌枕沉降量變化趨勢進行指數函數擬合，指數函數階數為2，擬合結果如圖11所示，則36 Hz激振力影響下軌枕沉降量與時間的演化規律大致服從圖中擬合過后得到的曲線發展規律。

分析處理得到道床累計下沉量隨時間變化的函數模型方程為

f(x)=-106×e-0.006 388x+106×e-0.000 638 7x

(1)

式中選取的常數項為最佳數值，即擬合度最高的沉降量計算式。當常數項在符合要求的區間時，可對相應的激振頻率進行擬合，可靠率達95%。

圖11 擬合曲線與仿真曲線對比圖

4 結束語

文中通過有限元軟件建立了縮小尺寸軌枕模型，得到網格劃分后的節點坐標，在離散元軟件中通過顆粒替換將軌枕單元顆粒按照節點坐標排列，得到有砟軌道結構的三維離散元模型。建立相應道砟箱模型，分析在穩定載荷下軌枕對道砟顆粒的作用情況、軌枕的受力情況和應力分布情況。

本研究根據仿真模型搭建了試驗平臺，綜合分析離散單元法在仿真道床實際狀態的優勢和特點，通過將模擬仿真分析結果與實際試驗數據對比，驗證了此耦合方法在穩定作業過程中仿真道床實際狀態的可行性和有效性，并利用離散單元法進一步深入研究。本研究發現軌枕在鋼軌與道床之間，上部穩定車傳遞過來的激振效果通過軌枕來影響道床，且與下部結構散體道砟顆粒直接接觸。軌枕正下方道砟顆粒振動速度均隨著道床深度的增加呈扇形迅速衰減，且隨著穩定作業的進行，軌枕下方道砟顆粒不斷向下移動，當激振頻率為36 Hz時，道床沉降量最大。利用EDEM軟件進行了作業過程的模擬仿真并提取了相關數據，得到了不同作業工況下任意時刻道床下沉量等相關數值，并給出了沉量隨時間變化的函數模型方程。