基于Hertz接觸理論的法向接觸剛度計算方法

關慶華 ，趙 鑫 ，溫澤峰 ，金學松

（西南交通大學牽引動力國家重點實驗室，四川 成都 610031）

車輛-軌道耦合動力學將車輛和軌道系統(tǒng)通過輪軌接觸來實現(xiàn)車輛系統(tǒng)和軌道系統(tǒng)的動態(tài)相互作用計算，真實地反映了車輛和軌道的動態(tài)特性，已廣泛應用于鐵道車輛和軌道的動力學計算、結(jié)構設計及安全評估等[1].其中，輪軌法向力的計算是確定車輛及軌道系統(tǒng)運動和振動的基礎，法向力的計算還將影響到接觸斑和輪軌蠕滑力的大小，進而將影響整個車輛和軌道系統(tǒng)的動態(tài)特性[2-3].

此外，通過有限元模型來分析車輛和軌道系統(tǒng)振動[4]以及計算高頻輪軌沖擊力[5]和輪軌噪聲[6]中，也需要精確的輪軌力確定方法.20 世紀70 年代，英國鐵路技術研究所的Jenkins 等[7]在研究車輛和軌道參數(shù)對輪軌垂向力的影響時，給出了錐形踏面和磨耗型踏面的輪軌接觸常數(shù)公式，分別如式（1a）和式（1b）.

式中：Gc和Gw分別為1∶20 錐形踏面和磨耗型踏面的接觸常數(shù)；參數(shù)R與車輪輪徑r相等.

該經(jīng)驗公式適用的輪徑范圍為0.150～0.600 m，對應的鋼軌頂面廓形的橫向曲率半徑為0.230 m.根據(jù)式（1a）和式（1b）可確定出接觸變形量δ和法向力P之間的關系為

孫翔[8]根據(jù)Hertz 彈性體接觸理論，針對輪軌接觸特定情況，給出了確定輪軌接觸橢圓的直接方法，可通過參數(shù)ρ/r（ρ為由輪軌主曲率之和確定的常數(shù)）直接確定輪軌接觸橢圓的各個參數(shù)及無量綱蠕滑系數(shù)，但未給出彈性接觸力及接觸剛度的計算方法.

彈性理論及接觸力學著作[2,9]中僅給出了Hertz彈性接觸斑長、短軸系數(shù)的計算公式，關于彈性體Hertz 接觸剛度的計算及其數(shù)表尚不多見，僅文獻[3]中給出了9 組插值數(shù)據(jù)，不便于工程計算.Hertz 接觸理論屬于經(jīng)典的彈性接觸力學，在實際的接觸過程中材料塑性及阻尼也會對接觸力產(chǎn)生影響[9-16].

本文基于Hertz 彈性體接觸理論[2-3]，給出了滿足Hertz 接觸條件的的彈性體接觸剛度通用計算公式，并針對輪軌接觸特定條件，在文獻[8]研究結(jié)果的基礎上，確定了輪軌接觸剛度的簡化計算方法和計算數(shù)表，以期對車輛-軌道耦合動力學中的輪軌力計算進行修正，并為輪軌摩擦磨損、高頻沖擊振動和噪聲研究中輪軌載荷確定提供計算依據(jù).

1 Hertz 彈性體接觸的基本公式

Hertz 接觸理論的基本假設如下：

條件1兩接觸體在初始接觸點附近的表面至少二階連續(xù).故表面在初始接觸點的二階微分鄰域內(nèi)可表示為橢圓拋物面、柱面或雙曲拋物面，可以用微分幾何的方法進行分析.

條件2接觸是非共形的.在加載過程中由初始接觸點逐漸形成接觸區(qū).由于接觸區(qū)很小，在初始接觸點附近兩彈性體可視為彈性半空間.

條件3小變形.在初始接觸點附近，兩接觸表面與初始接觸點處公切面法線相交的點相接觸.

條件4無面內(nèi)摩擦，因而切向面力為0.

一般而言，對于非共形輪軌接觸，條件1～3 是滿足的，但輪軌蠕滑的存在使得條件4 不滿足.Johnson[9]指出，兩非共形體接觸時，只有當其材料彈性常數(shù)不同時，摩擦力才會對法向接觸產(chǎn)生影響.因此，在目前的動力學分析中，通常認為法向接觸力和切向蠕滑力是解耦的，不考慮輪軌蠕滑對法向力的影響.

1.1 計算接觸橢圓及載荷變形關系的基本公式

考慮一般性條件，Hertz 接觸計算系數(shù)A和B可用式（3）、（4）求得.

式中：R11、R12和R21、R22分別為接觸體1 和接觸體2沿接觸坐標系縱向、橫向的主曲率半徑；ψ為兩接觸體包含R11和R21法平面的夾角.

根據(jù)Hertz 彈性體接觸理論，系數(shù)A和B與接觸橢圓偏心率e之間的關系為

式中：e2=1-（be/ae）2，ae、be分別為接觸橢圓長、短半軸長度；K（e）、E（e）分別為以Φ為積分變量的第一、二類完整橢圓積分：

接觸點確定后，即可通過接觸體主曲率半徑及夾角ψ，根據(jù)式（3）～（7）計算出接觸橢圓的偏心率及其對應的第一、二類完整橢圓積分.

根據(jù)接觸橢圓偏心率及第一、二類完整橢圓積分可確定出接觸橢圓長、短半軸長度分別為

式中：P為法向壓力；

其中：νi、Ei分別為接觸體i（i=1，2）的泊松比和彈性模量.

彈性接觸引起的壓縮量與法向力之間的關系為

式中：Knr為Hertz 接觸常數(shù)，可表示為

若換算為類似于式（1a）和式（1b）的表示方式，為

顯然，Knr和Gnr與彈性體的材料參數(shù)、主曲率半徑及其對應的橢圓完整積分有關.

1.2 以ρ/R11 為參數(shù)的接觸橢圓及載荷與變形計算

在輪軌接觸特定條件下，確定輪軌非線性接觸常數(shù)可進一步簡化.

鋼軌縱向曲率為0，即1/R21=0，同時忽略ψ后，式（3）和（4）分別為

式中：R11為車輪的法向滾動圓半徑，R11=rsecα，α為車輪踏面斜率角；R12為車輪踏面橫向曲率半徑，圓心指向車輪內(nèi)部時為正，否則為負；R22為鋼軌頂面橫向曲率半徑，凸形軌頭為正.

引入?yún)?shù)ρ：

式（18）給出了ρ/R11和θ的關系，分別對應接觸橢圓長軸垂直于軌道縱向和沿軌道縱向的情況.

由式（8）、（9）及式（14），可將ae和be以ρ/R11分別表示為

根據(jù)式（5）即可計算接觸橢圓的偏心率，對于輪軌接觸，可根據(jù)式（19）～（22），通過ρ/R11來確定接觸橢圓長短軸大小及Hertz 接觸常數(shù)Knr和Gnr.

為便于計算中查表使用，將θ、A/B、ρ/R11、m、n及qk制作成數(shù)表，如附加材料表S1 所示.其中，θ、A/B、m、n及qk為通用參數(shù)，適用于任何滿足Hertz接觸條件的彈性體接觸計算，而ρ/R11可用于輪軌接觸的查表計算，與文獻[8]的插值表格對應，可直接供車輛軌道耦合動力學計算使用.

需要指出的是：式（1a）和式（1b）中R為車輪半徑，即R=r，在文獻[1,7-8]中均以此計算.由于Hertz接觸是以接觸體的主曲率半徑來計算的，應以車輪的法向滾動圓半徑R11來代替，因此，本文的公式和數(shù)表均以ρ/R11為參數(shù)給出.

2 Hertz 非線性接觸剛度

Hertz 接觸常數(shù)Knr有時也被稱為Hertz 非線性剛度，單位為N/m2/3，因此，其并非傳統(tǒng)意義的剛度，它依賴于接觸體的材料參數(shù)、接觸點的主曲率半徑及兩接觸體包含R11和R21法平面的夾角，與載荷無關.

Jenkins 等[7]在計算輪軌沖擊力P1時，使用了線性化輪軌接觸剛度，其定義為

式中：khl與法向力P和靜輪載Po有關，需迭代求解.

根據(jù)剛度的定義，由式（10）可直接得到Hertz非線性接觸剛度為

顯然，khn不僅由Knr決定，還與P有關.

3 典型輪軌接觸分析

3.1 接觸參數(shù)qk 的計算分析

按照本文公式，計算得到qk與θ的曲線如圖1所示.圖中數(shù)據(jù)表明，θ越小，qk越大，對應的接觸常數(shù)Knr也越大；當θ＞30°時，qk隨θ的增大逐漸趨于定值.相比于文獻[3]中結(jié)果，本文數(shù)表中的qk覆蓋范圍更廣，利用本文公式可計算θ為任意值時的接觸常數(shù)Knr.

圖1 qk 與θ 的關系曲線Fig.1 Relationship curve of qk with respect to θ

3.2 典型輪軌接觸工況分析

圖2（a）所示為LM 車輪踏面及其橫截面曲率半徑.圖2（b）為CN60 鋼軌廓形及橫截面曲率半徑，軌頂圓弧半徑為0.300 m，向兩側(cè)依次為0.080 m和0.013 m 的圓弧.

圖2 車輪和鋼軌踏面Fig.2 Profiles of wheel and rail tread

以上述踏面為例，選取某地鐵輪軌參數(shù)：車輪名義半徑為0.365 m，輪緣內(nèi)側(cè)距為1.353 m，軌距為1.435 m，軌底坡為1∶40，對橫移量在-20～20 mm時的Hertz 接觸常數(shù)進行了計算，結(jié)果如圖3 所示.

圖3 LM 車輪踏面和CN60 鋼軌匹配下的接觸常數(shù)Fig.3 Contact constants for LM wheel profile and CN60 rail profile

由圖3 可知，應用本文式（11）和式（12）計算得到的Hertz 接觸常數(shù)Knr和Gnr相比于采用式（1a）和式（1b）兩種踏面計算公式的計算結(jié)果在橫移量-20～20 mm 內(nèi)均有較大差異.其中，基于式（1a）計算的Knr相比于本文結(jié)果的變化范圍為 -25.97%～71.23%，最小相差為0.68%；基于式（1b）計算得到的Knr均高于本文結(jié)果，相對變化范圍為0.40%～131.42%.基于式（1a）計算得到的Gnr相比于本文結(jié)果的變化范圍為 -30.13%～22.20%，最小相差0.45%；基于式（1b）計算的Gnr均低于本文結(jié)果，相對變化范圍為 -42.84%～-0.26%.

以名義接觸位置為例，利用式（11）計算的Knr為6.77 × 1010N/m3/2，而通過式（1a）和式（1b）計算得到的Knr分別為8.17 × 1010N/m3/2和11.08 ×1010N/m3/2，相對于式（11）計算結(jié)果分別增大了20.77%和63.77%.利用式（12）以及式（1a）和式（1b）計算的Gnr分別為6.02 × 10-8、5.31 × 10-8m/N2/3和4.33 × 10-8m/N2/3，利用式（1a）和式（1b）的Gnr計算結(jié)果相對于式（12）計算結(jié)果分別減小了11.82%和28.03%.

3.3 輪徑對接觸常數(shù)的影響

為分析輪徑對接觸常數(shù)Knr和Gnr的影響，分別對名義滾動圓位置和輪緣角接觸位置進行了計算.

圖4 為名義滾動圓附近計算結(jié)果，計算中，不考慮接觸角的影響，R12=-0.500 m，R22=0.300 m.圖中結(jié)果表明，名義滾動圓附近，Knr和Gnr精確計算結(jié)果與式（1b）計算結(jié)果相對接近，Knr的誤差范圍為0.02%～3.49%，與式（1a）計算結(jié)果相差較大，Knr的誤差范圍為25.94%～27.26%.

圖4 名義滾動位置的接觸常數(shù)Fig.4 Contact constants for nominal contact position

圖5 為輪緣角與鋼軌R13 圓弧接觸時的計算結(jié)果，計算中，接觸角取為70°，R12=-0.014 m，R22=0.013 m.圖中結(jié)果表明，輪緣角接觸位置，Knr和Gnr精確計算結(jié)果與式（1a）和式（1b）的計算結(jié)果相差均較大，其中，式（1b）計算Knr結(jié)果相對于本文式（11）計算結(jié)果的誤差范圍為17.32%～27.70%，式（1a）的計算誤差范圍為41.74%～44.89%.

圖5 輪緣角位置的接觸常數(shù)Fig.5 Contact constants for flange corner position

4 結(jié) 論

1）基于Hertz 彈性體接觸理論，給出了彈性體接觸參數(shù)的計算方法和接觸剛度參數(shù)數(shù)表，計算公式考慮到接觸體主曲率半徑、接觸斑形狀以及材料參數(shù)及載荷，計算數(shù)表可用于任意滿足Hertz 接觸條件（如輪軌、軸承、齒輪等）的彈性體接觸斑大小和接觸剛度的統(tǒng)一計算，填補了目前彈性理論及接觸力學中Hertz 接觸數(shù)表缺乏接觸剛度查表參數(shù)的不足.本文數(shù)表還融合了以ρ/R11為參數(shù)的輪軌特定接觸條件的接觸剛度插值，可為車輛-軌道耦合動力學中輪軌力計算提供修正.

2）目前采用的輪軌接觸常數(shù)計算公式為近似的經(jīng)驗公式，僅以錐形踏面和磨耗型踏面區(qū)分來考慮輪徑的影響，限定于特定的輪徑范圍和鋼軌廓形，在實際應用中具有明顯不足.以LM 踏面和CN60鋼軌的典型匹配為例，當車輪踏面名義中心附近圓弧與軌頂中心圓弧接觸時，現(xiàn)有經(jīng)驗公式中的磨耗型踏面公式計算結(jié)果與本文計算結(jié)果相比誤差較小；其它接觸位置時，經(jīng)驗公式與本文計算結(jié)果誤差較大，最大誤差超過131%.

3）本文的接觸剛度參數(shù)計算公式，適用于滿足Hertz 接觸條件的彈性接觸力計算.實際的接觸過程涉及復雜的力學問題和非線性因素，材料塑性變形、第三介質(zhì)、接觸能量耗散等都可能會影響到法向接觸力.在多體系統(tǒng)的接觸力計算中，除Hertz 彈性接觸力外，通常還引入阻尼項來考慮接觸過程中的能量損耗.此外，接觸過程中的塑性變形和黏彈性阻尼特性以及彈性體波動效應等也需要予以考慮.對于高速輪軌系統(tǒng)，如何考慮高速滾動接觸過程中的材料塑性、表面粗糙度、摩擦、彈性波動等引起的能量耗散對法向接觸力和蠕滑力的影響，尚有待于進一步的理論研究和試驗驗證.

備注：附加材料在中國知網(wǎng)本文的詳情頁中獲取.

致謝：感謝牽引動力國家重點實驗室自主研究課題（2020TPL-T02）支持.