賈昕昱,黎國清,劉秀波,楊 飛
(1.中國鐵道科學研究院,北京 100081;2.中國鐵道科學研究院集團有限公司 基礎設施檢測研究所,北京 100081)
滾動接觸疲勞裂紋是重載線路上最為常見的一種傷損形式,這種裂紋的萌生和擴展與輪軌之間的接觸狀態有很大關系。李霞[1]分析了鋼軌軌底坡與LM,LMA型2 種車輪的接觸狀態。Vo K D[2]分析了鋼軌存在磨耗時的輪軌接觸狀態。但是對于輪軌接觸狀態與裂紋萌生之間關系的研究還相對較少。
有限元方法是分析輪軌接觸狀態非常有效的一種方法,三維瞬態滾動有限元模型可以得到車輪、鋼軌的受力狀態以及輪軌接觸斑的應力分布等信息。于淼[3]采用三維瞬態有限元模型分析了車輛經過波磨區段時輪軌間的動力響應特性。在對裂紋萌生的分析中,周宇[4]采用多體動力學軟件得到輪軌接觸斑的應力分布狀態,再將接觸斑內的力施加到鋼軌上,用以對裂紋萌生狀態進行分析。但是這些分析過程忽略了車輪滾動過程中接觸斑的分布狀態是不斷變化的。
本文基于臨界平面理論建立鋼軌疲勞裂紋萌生預測模型,采用三維瞬態滾動有限元模型模擬重載鐵路的輪軌接觸狀態,并分析不同軌底坡和輪軌摩擦系數下輪軌間的接觸狀態及鋼軌裂紋萌生特征,得到不同因素對裂紋萌生的影響。
對于材料的疲勞裂紋萌生壽命預測,目前大體可分為2 種方法,一種是通過實際測量和模擬試驗得到疲勞數據,但是這種方法沒有較好的通用性,另一種是疲勞壽命分析法,通過考慮材料的載荷形式以及材料的相應參數分析結構的疲勞壽命,而臨界平面法就是這種方法中較為有優勢的方法。臨界平面法是根據剪應力以及拉應力確定一個物理量,找到這個物理量的最大值所在的平面即臨界平面,材料在臨界平面上最先產生裂紋。
Jiang 和Sehitoglu[5-6]提出的疲勞參量FP,是基于應變能與臨界平面法的一個物理量,為

其中,

在有限元分析中,OXYZ為全局坐標系;O′X′Y′Z′為旋轉平面上的局部坐標系,并且分別表示鋼軌的縱向、橫向和垂向,θ為裂紋所在平面的法向向量N與Z′軸間的夾角,?為裂紋所在旋轉平面的法向向量在X′Y′平面內的投影N′XY與X′軸的夾角,如圖1所示。

圖1 旋轉平面示意圖
通過有限元分析可以得到所求節點的應力分量σij和應變分量εij,根據張量分析理論可以組成相應的應力矩陣σ和應變矩陣ε。將所得的應力、應變矩陣進行旋轉變換,就可得到任意角度平面上的應力分量σ′ij和應變分量ε′ij,分別為

其中,

式中:A為張量旋轉變換矩陣。
通過對相應單元的疲勞參量進行分析,找到疲勞參量FP的最大值FPmax,其所在平面就是疲勞裂紋萌生的平面,也就是臨界平面。這一模型可以將拉應力、剪應力、拉應變、剪應變有機地結合起來,充分考慮到拉應力和剪應力同時對裂紋萌生和擴展的影響。
Akama[7]提出拉伸型裂紋壽命預測模型為

式中:σ′f為拉伸強度系數;E為彈性模量;Nf為輪對的通過次數;ε′f為拉伸延性系數;b為疲勞強度系數;c為疲勞延性系數。
王建西[8]采用剪切型裂紋壽命預測模型,為

式中:τ′f為剪切強度系數;G為切變模量;γ′f為剪切延性系數。
在分析車輛與軌道的相互關系時,傳統的多體動力學軟件往往傾向于分析結構相對低頻的響應,對于高頻響應有一定的缺失。而鋼軌裂紋的產生和擴展更傾向于輪軌之間的高頻作用,傳統的多體動力學軟件無法較好地解決高頻振動問題,因此,采用ANSYS/LS-DYNA 軟件建立輪軌三維瞬態滾動接觸有限元模型,通過該模型可以分析重載列車與有砟軌道之間的相互關系。由彈簧質量、一系懸掛、輪對、鋼軌、扣件、軌枕和道砟組成的有限元模型及輪軌尺寸如圖2所示。
由于一系懸掛以上的車輛部分受到高頻振動的影響較小,因此車輛一系懸掛以上的部分被簡化為1 個剛性質量塊,質量塊的大小根據列車的軸重來決定。一系懸掛由8 根并列的彈簧阻尼結構組成,用以連接質量塊及車軸。由于輪對滾過的距離較短,仿真時忽略滾動過程中輪對橫移的影響,但對輪對的橫向運動進行約束。車輪及軌道采用真實重載線路中的參數進行模擬,車輪踏面為LM 磨耗型,鋼軌為75 kg·m-1的U75V 熱處理鋼軌。車輪與鋼軌之間的接觸采用基于罰函數的面—面接觸準則,接觸面之間的摩擦采用庫侖準則。

圖2 輪軌三維瞬態滾動接觸有限元模型
圖2(b)中,A 點為輪對滾動的起始點,在該點計算輪軌系統的初始變形。從A點到B點為輪對的動態松弛區,以保證輪對能夠平穩地滾動,距離為0.35 m。從B 點到C 點為疲勞裂紋萌生計算區,用以計算鋼軌疲勞裂紋的萌生。軌道的總長度為18 m,其間包括30 個軌枕。在動態松弛區以及計算區內鋼軌單元劃分的較為細密,遠離輪對滾動的地方單元劃分較為稀疏。為了得到較為精確的輪軌接觸狀態,在車輪與鋼軌可能的接觸區域內均劃分較為細密的單元,在輪對的起始滾動點即A 點處,鋼軌的最小網格尺寸為1 mm,在動態松弛區內最小網格尺寸為2.5 mm,在疲勞裂紋萌生計算區內,最小網格尺寸為0.7 mm。整個模型的總節點數和總單元數分別為1.18×106個和1.04×106個。有限元模型中所用參數見表1[9]。
有限元分析得到車輪滾過平順鋼軌時的輪軌垂向力如圖3所示。圖中:A 點,B 點和C 點為圖2(b)中的位置點。從圖3可以看出:從A 點到B 點的范圍內,輪軌垂向力的波動幅值較大,沒有達到穩定的狀態,這就是由于輪對從靜態到動態轉變時不可避免地會產生振動;隨著時間的增加,在阻尼的作用下,輪對的振動幅值在逐漸減小,當到達B點時,輪對基本處于穩態滾動狀態;根據Zhao 等人[10]研究的結論,當輪軌垂向力的波動幅值除以輪對的靜軸重小于10%時,即可認為車輛達到穩態運動狀態,有限元結果即可用于相應的計算分析,從B點到C點,輪軌垂向力基本平穩,其最大波動幅值除以單個車輪的靜軸重約為5.5%,且輪軌垂向力在125 kN 附近波動,與軸重相吻合,表明所建模型可以用來進行進一步的分析計算。

表1 有限元模型中所用參數

圖3 輪軌垂向力
由于不同軌底坡會導致輪軌接觸點的位置及接觸斑幾何形態的改變,這一改變可能導致疲勞裂紋萌生位置及萌生時間的改變,因此,分析不同軌底坡時的輪軌接觸狀態對預測疲勞裂紋的萌生起著重要的作用。我國的重載鐵路中,軌底坡多數設置為1∶40,但是在一些早期的線路中,仍存在軌底坡為1∶20 的情況,因此考慮軌底坡分別為1∶40,1∶30,1∶20和1∶10的情況。
首先,分析輪對相對于軌道發生一定范圍橫移時的輪軌接觸狀態。輪背內側距為1 353 mm,軌距為1 435 mm,以此確定車輪與鋼軌的初始接觸位置。以輪軌初始接觸位置為起點,輪對分別向左右移動,移動的最大距離為12 mm,并用直線將車輪上的點與鋼軌上的對應點相連,軟件截屏得到整個移動范圍內的最短接觸長度如圖4所示。圖中:紅線為輪對沒有橫移時的輪軌接觸狀態。

圖4 不同軌底坡時輪軌接觸狀態(單位:mm)
從圖4可以看出:軌底坡為1∶40 時,輪軌接觸位置主要位于鋼軌的軌距角處,集中在距軌距邊9~22 mm 范圍內;軌底坡為1∶30 時,輪軌接觸位置主要還是集中在軌距角處,但是在軌頂處也出現了少量的分布,且集中在距軌距邊12~27 mm范圍內;軌底坡為1∶20 時,輪軌接觸位置在軌頂處及軌距角處均有分布,并且主要集中在軌頂的兩邊,距軌頂2 側13 mm 范圍內;軌底坡為1∶10時,輪軌接觸位置在軌頂外側,距軌距邊16~21 mm范圍內。可以看出,改變軌底坡將改變輪軌接觸位置,而軌底坡為1∶40 時輪軌接觸位置主要集中在軌距角處,這一情況也有可能是鋼軌斜裂紋的萌生因素,需要做進一步的分析。
通過分析車輪與鋼軌在初始接觸位置沿著軌道直線運動的情況,得到不同軌底坡時鋼軌的壓力分布如圖5所示,從而得到不同軌底坡時的法向節點力和接觸斑形態,如圖6所示。圖6中:上半部分為接觸斑內的法向節點力,映射到下平面中的為接觸斑的形態。

圖5 不同軌底坡時鋼軌壓應力分布
從圖5和圖6可以看出:軌底坡為1∶40 時,輪軌接觸斑不在鋼軌軌頂處,而是偏向于軌距角處,接觸斑呈橢圓形,其長半軸和短半軸的長度分別為9.0 和4.5 mm,最大壓應力為1.20 GPa;軌底坡為1∶30 時,接觸斑的位置比軌底坡為1∶40時更靠近軌頂,接觸斑的面積增大,不再是標準的橢圓形,最大壓應力略微降低,為1.16 GPa;軌底坡為1∶20 時,接觸斑不再是1 個橢圓形,而變成了1 個帶狀分布的圖形,整個接觸斑位于軌頂處,在軌頂2 側有2 個峰值,最大壓應力大幅降低,為0.70 GPa;軌底坡為1∶30時,接觸斑從鋼軌軌頂內側轉移到鋼軌軌頂外側,最大壓應力為1.31 GPa;當軌底坡為1∶20時,車輪和鋼軌發生了共形接觸,導致接觸面積增大,同時壓應力大幅降低,這種軌道參數與LM型車輪更為匹配。

圖6 不同軌底坡時法向節點力和接觸斑形態
最后,為了驗證上述的輪軌接觸狀態,仿真得出不同軌底坡時車輪與鋼軌的接觸區域,且由于輪軌接觸區內單元的變形為0.05 mm,僅選擇車輪與鋼軌接觸距離小于0.05 mm 的線段,結果如圖7所示。從圖7可以看出:圖(a),圖(b)和圖(d)中的輪軌接觸區域均靠近軌距角處;圖(c)中均勻分布在軌頂的2 側,且接觸區域遠大于圖(a),圖(b)和(d),這一結果也與圖5中的結果相一致。

圖7 不同軌底坡時輪軌接觸區域
為了找到疲勞裂紋萌生平面,將應力應變矩陣進行旋轉,旋轉矩陣A中θ和?的變化范圍均為0°~360°,步長為5°。分析軌底坡為1∶40、摩擦系數為0.5 的情況,提取出計算區域內各單元隨時間變化的應力應變數據,代入式(1)中進行計算,求出各單元在不同旋轉角度下的疲勞參量如圖8所示。

圖8 各單元疲勞參量
圖8中左側縱坐標數值為式(1)中的拉應力/拉應變部分;右側縱坐標數值為式(1)中的剪應力/剪應變部分。從圖8可以看出:當疲勞參量FP小于0.1 時,有部分單元是由拉應力及剪應力共同作用;當FP大于0.1 時,單元受力以剪應力為主;圖8中,當單元受到壓應力的作用時,式(1)中拉應力/拉應變部分的值為0,疲勞參量值由剪應力/剪應變部分決定;當單元受到拉應力的作用時,式(1)中剪應力/剪應變部分近似為0,疲勞參量值由拉應力/拉應變部分決定。在各單元中,疲勞參量的最大值為0.67,拉應力/拉應變為0,因此壽命預測公式選擇剪切型裂紋壽命預測模型,即式(7)。
在單元疲勞參量為最大值時,對應的旋轉角度為θ=40°,?=10°。將所得角度與在實際線路中觀察到的裂紋形態進行對比,結果如圖9所示。圖中:行車方向水平向右,裂紋面的垂線用虛線表示,與X軸的角度為θ,與Z軸的角度為?。從圖9可以看出,計算結果中的角度與實際裂紋的形態較為一致。

圖9 裂紋角度
分析得到各單元的疲勞參量如圖10所示。圖中:紅線輪廓為軌頭廓形,柱狀圖的長度為疲勞參量。從圖10可以看出:疲勞參量最大值發生于軌頂以下1.96 mm 處,到軌距角的水平距離為21.5 mm,這一結果與現場觀測到的斜裂紋基本一致。

圖10 各單元的疲勞參量
對不同輪對橫移量時裂紋的萌生進行分析,將輪對向外軌方向進行橫移,橫移量分別為6 和12 mm,輪對橫移后的輪軌接觸狀態如圖11所示。從圖11可以看出:當輪對橫移量為12 mm 時,接觸斑位于軌距角處,接觸班呈長條狀,最大壓力為1.83 GPa;輪對橫移量為6 mm 時,接觸斑呈近似圓形,最大接觸壓力為953 MPa。

圖11 不同輪對橫移量時輪軌接觸狀態
不同輪對橫移量時對應的疲勞參量如圖12所示。從圖12可以看出:當輪對向外軌側橫移時,疲勞參量在增大,更容易導致裂紋的萌生。

圖12 不同橫移量的疲勞參量
由于疲勞參量為最大值時,對應單元的拉應力/拉應變部分為零,只存在剪應力/剪應變部分,因此采用剪切型壽命預測公式,即式(7)進行裂紋萌生壽命分析。材料參數參考Jiang[11]所提供的數據,式中的τ′f,γ′f,b和c分別取為1 015,0.509 4,-0.104 8 和-0.550 1。計算得出在這種情況下,鋼軌疲勞裂紋萌生時輪對的通過次數Nf為9.19×105次,通過總重為2 298 萬t。由于實際輪軌接觸位置是在1 個接觸區域內隨機變化的,而文中所分析的接觸區域是相對固定的,因此實際的萌生壽命比文中預測的要長,文中得出的計算結果是偏于安全的。
采用相同的方法,計算出軌底坡為1∶30,1∶20 和1∶10 時鋼軌疲勞裂紋的萌生壽命分別為3 934 萬t、9.88 億t 及3 875 萬t,如圖13所示。從圖13可以看出:軌底坡由1∶40 增加到1∶20 時,鋼軌疲勞裂紋的萌生壽命增加;軌底坡由1∶20 增加到1∶10 時,萌生壽命降低,而且軌底坡為1∶20 時的萌生壽命最大,說明LM 型車輪在軌底坡為1∶20 時與75 kg·m-1鋼軌更為匹配,此時可以大幅提升裂紋的萌生壽命,提高鋼軌的使用時間。

圖13 不同軌底坡裂紋萌生壽命
1∶40 軌底坡時不同摩擦系數對鋼軌疲勞裂紋萌生壽命的影響結果如圖14所示。從圖14可以看出:隨著摩擦系數的增大,通過總重隨之減小,也就是摩擦系數越大越容易導致裂紋提早萌生。

圖14 不同摩擦系數對鋼軌疲勞裂紋萌生壽命的影響
為了探究這一原因,分析不同摩擦系數下接觸斑內節點的最大法向力和最大切向力,結果如圖15所示。從圖15可以看出:節點最大法向力隨摩擦系數的變化不大,節點最大切向力隨摩擦系數增大而增大,也就是法向力對裂紋萌生不起主要作用,而切向力越大越容易導致裂紋盡早萌生,切向力是裂紋萌生的主要影響因素。

圖15 不同摩擦系數時節點最大法向力和最大切向力
(1)軌底坡分別為1∶40,1∶30,1∶20 和1∶10 時,隨著軌底坡的增大,接觸斑的位置先逐漸向軌頂發展,再向軌頂外側發展,接觸斑的面積先逐漸增大,再降低。車輪踏面為LM 型時,在1∶20軌底坡時與75 kg·m-1鋼軌更為匹配。
(2)軌底坡為1∶40 時,疲勞裂紋的萌生位置位于軌頂下方1.96 mm 處,到軌距角的水平距離為21.5 mm,這一結果與現場觀測到的斜裂紋基本一致,證明了文中所使用的模型可以有效預測鋼軌斜裂紋的萌生。
(3)軌底坡分別為1∶40,1∶30,1∶20 和1∶10 時,隨著軌底坡增大,鋼軌疲勞裂紋萌生時的通過總重先增大再降低,在軌底坡為1∶20 最大,1∶20的軌底坡可以減緩疲勞裂紋的萌生。
(4)摩擦系數分別為0.25,0.30,0.35,0.40 和0.50 時,隨著摩擦系數增大,節點法向力基本不變,節點切向力逐漸增大,疲勞裂紋的萌生壽命逐漸降低,表明摩擦系數越大,切向力越大,疲勞裂紋越容易萌生。