飛機輪組滑水行為與道面積水附加阻力

李岳，宗輝杭，蔡靖，戴軒

(中國民航大學 交通科學與工程學院，天津 300300)

根據國際航空運輸協會統計，近十年全球共發生708 起重大航空安全事故，飛機沖偏出跑道事故占比達23%，在全部事故類別中排名第一[1]，有半數此類事故伴隨道面積水情況。飛機地面高速滑行時道面積水可引起輪胎附著力與附加阻力降低，誘發單個輪胎滑水乃至整機滑行姿態失控。因此，飛機輪組滑水行為與附加阻力是滑水問題研究的重點。

美國國家航空航天局（NASA）較早開展了飛機輪胎滑水性能試驗，在此基礎上建立臨界滑水速度（后文簡稱vP）計算公式[2]（NASA 公式）。一般認為輪胎接觸面上動水壓力的豎向抬升作用是導致滑水的主要原因，后續研究中廣泛采用道面對輪胎豎向支撐力作為輪胎滑水判定指標[3-8]。近年來，有限元軟件和流固耦合算法的成熟使得仿真分析逐步成為滑水研究的常見手段。Ong 和Fwa[3-4]通過輪胎與積水路面流固耦合分析，驗證當流體對輪胎豎向支撐力與車輪荷載平衡時即發生打滑現象。Fwa 和Ong[5]研究了胎面下部流體抬升力對輪胎接觸面法向接觸力的削弱效應。Tang 等[6]基于實際路面掃描數據建立仿真分析模型，探討單個輪胎在不同降雨強度和路面滲透性下的滑水現象。Srirangam 等[7]采用單輪仿真分析模型計算不同滑移率下的vP，分析路面積水引起的制動牽引力和轉向力損失。劉修宇等[8]開展輪胎滑水行為仿真，分析滑水前后輪胎路面豎向接觸力變化。黃曉明等[9]通過ABAQUS 仿真回歸得出了水流豎向托舉力與水膜厚度和行駛速度的關系式。鄭彬雙等[10]建立單輪滑水分析模型，研究濕滑狀態下瀝青路面附著系數，探討滑移率對輪胎滑水的影響。閆珅[11]提出了道面污染物阻力的計算公式并對污染跑道條件下的飛機起降性能進行分析。李岳等[12]開展輪胎與積水道面流固耦合分析，推導了輪胎接觸面動水壓強與道面豎向支撐力表達式。Fwa 等[13]通過有限元仿真探討臨界滑水速度與路面車轍深度之間的關系。

綜合來看，以往研究工作主要關注單輪滑水形成過程，將道面支撐力下降至0 時刻作為滑水臨界狀態，假定滑水過后輪胎與道面摩擦力完全消散，而對輪胎所受道面積水附加阻力討論比較少見。由于道面支撐力與飛機行進方向正交，而積水附加阻力與飛機行進方向一致或處于同一運動平面內，滑水臨界時刻積水附加阻力曲線出現反彎點，符合滑水狀態判定需要，因而積水附加阻力更適合作為滑水分析指標。當考慮瀝青道面不均勻積水條件時，輪組內部存在滑水過程差異，無法用單輪滑水行為表征。

本文以積水附加阻力為切入點，研究單個飛機輪胎與輪組滑水行為，分析積水附加阻力形成原理，建立流固耦合仿真模型，探討不同滑行狀態與積水條件下附加阻力規律，討論輪組附加阻力特征與疊加條件，為輪胎滑水條件判定和滑行姿態控制提供分析依據。

1 積水附加阻力構成

道面積水主要通過擠壓和噴流2 種方式對地面高速滑行滾動輪胎產生阻力效應，削弱輪胎與道面接觸條件，降低飛機起降性能，增大滑水事故發生幾率。英國布里斯托爾大學通過室內試驗研究了飛機輪胎濺水行為[14]，輪胎濺水噴流如圖1所示。

圖1 飛機輪胎濺水噴流試驗Fig.1 Water spread test of aircraft tire

輪胎駛過積水道面時積水與輪胎前緣接觸，因擠壓效應形成位移阻力分量Ddisp；部分流體向下通過底部胎紋構成的排水通道向后排出，輪胎接觸面上摩擦阻力分量Dfric隨速度增加而逐步減弱；部分流體在擠壓作用下橫向繞流，與輪胎邊緣阻擋水墻沖撞形成側方噴流；其中指向兩側的低強度噴流水柱較弱，對飛機影響有限；指向側后的高速噴流擊中機身或起落架時會產生沖擊阻力分量Dimpa；噴濺到機體上的水沿著表面流動繼而形成表面阻力分量Dsurf，即有道面積水總附加阻力為

Dimpa和Dsurf由濺水噴流引起，與機身外形及投影面積相關。文獻[15]的研究表明，前起落架是Dimpa和Dsurf的主要來源，主起落架濺水噴流方向相對遠離機身，當積水厚度小于6.5 mm 時濺水噴流阻力效應微弱[15]。輪胎發生滑水后，Dfric近似消散，因而可以將Ddisp作為輪胎附加阻力代表值，該指標是本文滑水分析的重點。主起落架輪胎位移阻力為[16]

式中：ρ為積水密度；vG為 地面滑行速度；n為輪組系數，典型雙輪起落架可近似取2.0，即為單輪式的2 倍。阻力系數CD可通過有限元仿真計算得出，輪胎發生滑水前近似取保守值0.75，當vG>vP時需考慮輪胎滑水影響。輪胎前緣積水接觸面面積S為[16]

式中：d為積水厚度；b為接觸面有效寬度；W為輪胎寬度；δ為輪胎垂直形變，可由輪胎“載荷-變形”曲線查得。輪胎發生滑水時，表層水膜將輪胎與道面完全分隔，圖1 中濺水噴流路徑隨之改變，滑水速度影響因子F為

由式（2）可知道面積水附加阻力隨滑行速度增加而增大，在vG/vP<1 階段，F保持為1.0，無修正作用，此時附加阻力僅為關于vG的二次遞增函數，vG/vP=1.0時達到最大值；滑水速度對附加阻力影響主要體現在vG/vP>1.0階 段，F隨滑行速度增加而出現非線性衰減；圖2 中滑水臨界狀態前后附加阻力變化規律有明顯轉折，可作為滑水判定依據，用于飛機滑水分析。

圖2 滑水速度影響因子變化曲線Fig.2 Variation curve of impact factor of taxing speed

2 飛機輪組滑水分析模型

選取空客A320 飛機為案例，采用ABAQUS 軟件建立主起落架雙輪流固耦合滑水分析模型。結合文獻[12]中的研究經驗，選取水流沖擊方式模擬輪組與積水道面相互作用，計算效率及精度可滿足分析需要[17]。

2.1 模型參數

基于實測數據繪制主起落架輪胎橫斷面，規格為46×17R20，繞輪軸旋轉一周生成三維實體模型，斷面細部尺寸如圖3 所示。

圖3 機輪模型生成過程Fig.3 Generation process of aircraft tire model

輪胎模型關于輪轂中心O點剛性耦合，并在該點施加輪載154.5 kN，通過調整輪載標準值考慮高速著陸條件下升力對飛機輪胎滑水的影響；輪壁內側施加1 140 kPa 均勻分布的壓力，模擬輪胎充氣效果；查詢飛機設計手冊雙輪中心間距0.93 m。輪胎結構等效為具有統一超彈性特性的均質橡膠材料，采用Mooney-Rivlin 本構模型模擬[18]，相關參數如表1 所示[19]。

表1 輪胎材料物理力學參數[19]Table 1 Mechanical parameters of tire material[19]

跑道接地帶長期承受機輪反復碾壓作用，容易發生輪轍變形（見圖4（a）），形成局部積水增厚和不均勻分布，對輪組阻力影響不容忽視。對此，參照美國國家道面測試中心試驗建立帶輪轍道面分析模型，對照平整道面情況。該試驗對瀝青道面進行飛機輪組足尺循環加載，模擬加載位置橫向偏移，輪轍斷面更符合實際道面變形特征[20]。

圖4 道面模型斷面特征Fig.4 Cross-sectional feature of pavement model

圖4（b）中輪轍斷面總影響寬度超過9.0 m，而飛機輪組寬度僅為136 cm，輪轍影響寬度達到飛機輪組寬度的6 倍以上，最大輪轍深度為12.1 mm，對此定義3 組加載荷位開展對比分析。其中，A1 加載荷位對應輪轍橫斷面最大深度處，輪組整體積水深度較大；A2 加載荷位對應輪轍橫斷面斜率較大處，左右輪胎之間積水深度差較大；A3 加載荷位對應輪轍橫斷面平坦段，輪組內部積水深度差較小，但整體積水厚度仍高于平整道面情況。3 組荷位涵蓋了跑道接地帶典型積水特征區域，可用于開展輪組滑水仿真分析（注：圖中橢圓形僅為輪組加載荷位示意，非成比例繪制）。

裝配完成后的輪組滑水仿真模型如圖5 所示，道面平面尺寸經反復試算確定為2.0 m×1.5 m，對圖4 中不同輪轍特征區分段進行模擬；根據跑道允許起降條件確定積水層厚度范圍為3～13 mm[12]，且在積水層上部定義300 mm 的厚空氣層，為濺水噴流提供充分的運動表達空間。積水層物理參數在文獻[21]中已詳細說明，此處不再贅述。

圖5 飛機輪組滑水仿真模型Fig.5 Simulation model of hydroplaning analysis of aircraft landing gears

采取耦合的歐拉拉格朗日（coupled Eulerian-Lagrangian，CEL）算法求解流固耦合問題，飛機輪組由8 節點線性六面體減縮積分單元（C3D8R）離散，道面采用剛性殼單元模擬，滑水過程中道面無變形發展；積水和空氣層采用歐拉網格離散（EC3D8R 單元）。仿真開始后水流由入口側按一定速度沖擊轉動輪胎形成附加阻力，在輪胎周圍產生濺水噴流，隨后由壓力出口側排出；調節水流速度使輪胎逐步達到臨界滑水狀態。受道面輪轍影響，各個輪胎運動協調關系相互獨立，保證其與道面充分接觸；輪胎中心轉動角速度與水流速度匹配，滑移率保持為1.0。

2.2 模型驗證

輪胎剛度特征多采用豎向位移量或輪印面積指標進行驗證，在文獻[19，21-24]中被廣泛采用，而輪胎滑水行為特征是模型驗證的難點。在ABAQUS分析中定義2 組歷程輸出，分別提取輪胎迎水面上位移阻力Ddisp和道面對輪胎豎向支撐力的結果，如圖6 所示。

圖6 輪胎模型滑水特性驗證Fig.6 Verification of hydroplaning feature of tire model

圖6 中支撐力曲線隨滑行速度增加逐步下降，支撐力降至0 時速度約為202 km/h，與NASA 公式計算該胎壓下vP=214 km/h 較為接近。位移阻力以散點形式繪出，受濺水噴流沖擊效應影響，Ddisp先振蕩增大后逐步減小，變化規律與由式（2）和式（4）繪制的經驗公式曲線相一致，且在208 km/h 附近達到最大值，與前2 組結果相差5%以下。因而，本文建立的仿真模型滑水特性得到驗證，可再現輪胎滑水前后位移阻力的變化過程。

值得注意的是，當滑行速度超過vP時道面支撐力保持為0，Dfric完全消散，而Ddisp仍處于非線性下降階段，持續影響飛機姿態和滑行距離，Ddisp更適合作為滑水狀態判定依據。

3 積水阻力結果分析

實際道面起降環境復雜，依托本文建立的飛機輪組滑水分析模型，開展不同飛機滑行狀態與道面積水分布條件下的積水阻力分析，探討輪組阻力特征與疊加關系。

3.1 飛機滑行狀態影響

飛機著陸過程中輪胎瞬間高速接地開始滾動滑行，在發動機反推和剎車系統共同作用下持續減速直至脫離跑道，輪胎受到積水阻力規律與一般起飛加速過程有明顯差異。有統計表明，著陸階段發生滑水事故幾率達到起飛階段的4 倍[25]，滑行狀態影響不容忽視。

假定道面覆蓋7.6 mm 厚均勻積水，水流以250 km/h 初始速度沖擊輪組模型，當輪胎接觸面水域完全覆蓋并穩定后，降低速度模擬著陸滑行過程，作為起飛加速情況對照，其他參數與2.1 節一致。圖7 為不同滑行速度下著陸滑行道面積水分布。

圖7 著陸滑行過程道面積水分布Fig.7 Distribution of accumulated water on pavement surface during landing

圖7 中上層Eulerian 網格描述流體運動分布，下部規整道面網格為流體變形分析提供參考系。如圖7（a）所示，著陸瞬時飛機滑行速度超過vP，輪胎底部可見大量水流通過，積水承托作用使得輪胎與道面脫離，輪胎接觸面輪廓不清晰，處于完全滑水狀態；當速度降至210 km/h 時，如2.2 節討論，此時輪胎處于滑水臨界狀態，輪胎接觸面恢復且較為完整，胎紋為流體提供底部排水通道，輪胎前緣及兩側Eulerian 網格“擠密”，積水阻力處于較高水平；隨著滑行速度進一步降低，圖7（c）中流域已被左右機輪接觸面完全分隔，輪胎脫離滑水狀態，除輪胎前緣出現雍水外，輪胎外側及中間輪隙Eulerian 網格變形明顯，繞流流體與阻擋水墻沖撞形成側方噴流。顯然，不同于飛機起飛階段由低速充分接觸發展至高速滑水脫離過程，著陸滑行狀態可改變輪胎與道面接觸條件及積水分布規律。

圖8 為不同滑行狀態下道面積水阻力，依據式（2）和式（4）對散點數據分段擬合，并以擬合曲線形式繪出。如圖中箭頭指向，著陸滑行時輪胎所受道面積水阻力經歷了反向先增大后減小的變化過程，阻力變化路徑與起飛滑行并不重合。相同速度條件下著陸滑行阻力更小，峰值點對應vP更低，超越臨界滑水狀態概率增加，符合滑水事故統計規律。積水厚度為7.6 mm 時，vP結果相差約8.3%，積水厚度為13.0 mm 時，vP相差增大至10.6%，與文獻[25]中給出的15%參考值相接近。

圖8 不同滑行狀態下道面積水阻力Fig.8 Drag force of accumulated water on pavement surface under different taxing conditions

綜合本節所述，飛機滑行狀態對輪胎滑水影響不容忽視，積水阻力結果表明：著陸時滑水風險更高，滑水控制條件更嚴格，后續分析針對著陸滑行條件展開。

3.2 道面積水條件影響

道面積水分布是影響輪胎滑水行為的另一重要因素。當局部輪轍變形形成后，可干擾道面橫坡正常排水功能，增大道面積水阻力，引發輪組內部滑水行為差異。

假定環境降雨與道面排水達到動態平衡，流域內平均積水厚度為7.6 mm。以A1 加載荷位為例（見圖4），此時輪轍中心最大水深可達19.7 mm，發生局部積水厚度超限。依次對圖4 中3 組加載荷位進行滑水分析，每組加載荷位具有獨立的輪轍變形與積水特征，并與平整道面情況對照。其他參數與2.1 節保持一致。機輪作用于不同荷位時輪組前緣動水壓強分布如圖9 所示。

圖9 輪胎前緣動水壓強分布Fig.9 Distribution of hydrodynamic pressure at tire frontier

由圖9 可知：①平整道面條件下積水均勻，輪胎前緣動水壓強關于輪組中心點對稱分布，左右輪胎壓強峰值同時達到胎壓等級1.14 MPa，積水阻力擠壓輪胎變形，可視為輪胎滑水發生臨界時刻，輪組內部滑水進程基本一致；②當輪組作用于A1 加載荷位時右輪積水厚度高于左側，右輪達到滑水臨界時刻時，左側壓強峰值僅為0.92 MPa，而左輪達到滑水臨界時刻時，右側壓強峰值已降為0.61 MPa，對照圖8 中滑水前后積水阻力變化規律，兩側輪胎并非同時發生滑水；③A2 加載荷位位于輪轍斷面大斜率段，整體積水厚度雖低于A1 加載荷位，但兩側積水厚度差加大，左輪積水更厚因而可先于右輪達到滑水臨界時刻，圖9 中動水壓強曲線非對稱分布特征顯著；④A3 加載荷位位于輪轍斷面相對平坦段，整體積水厚度和兩側厚度差均小于A1 和A2加載荷位，動水壓強分布非對稱性減弱，較平整道面仍有差異。

圖10 為不同加載荷位時左右側輪胎積水阻力變化規律，結果以擬合曲線形式繪出。表2 不同加載荷位臨界滑水速度比較結果?？梢钥闯觯孩僭?50 km/h以下飛機著陸滑行速度范圍內，3 組加載荷位道面積水阻力均高于平整道面情況，阻力曲線上升與下降速率增大，vP結果低于平整道面情況，對飛機滑行安全不利；②除平整道面情況外，3 組加載荷位輪組阻力曲線并不重合，兩側輪胎并非同時發生滑水，符合圖9 中動水壓強分布規律；③A1 加載荷位覆蓋輪轍斷面最大深度處，vP結果與平整道面相差最大達到12.8%，且右輪vP較左輪低3.5%；④與A1 加載荷位相反A2 加載荷位左輪vP低于右輪，符合輪轍斷面局部變形與積水特征，由于積水厚度差高于其他荷位，輪組內部vP相差達到5.7%；⑤A3加載荷位輪轍變形相對平坦，整體積水厚度略高于平整道面，因而輪組內部vP相差僅為1.1%，較平整道面下降6.8%。

表2 不同加載荷位臨界滑水速度比較Table 2 Comparison of hydroplaning speed at different loading positions

圖10 不同加載荷位時道面積水阻力Fig.10 Drag force of accumulated water on pavement surface at different loading positions

顯然，道面積水條件對輪組滑水影響不容忽視，其內部滑水過程存在時空差異，輪胎前緣動水壓強非對稱分布，積水阻力變化曲線相互獨立。輪轍引起局部積水增厚導致vP降低，滑水風險隨之加大。

3.3 輪組阻力特征

對于輪組滑水問題，一般認為由于主起落架輪胎物理間距較大，積水濺水噴流干擾效應有限，對于理想平整道面情況，輪組滑水可近似由單輪滑水行為代表，當起落架構型較簡單時，輪組阻力約等于單輪阻力的疊加，式（2）中輪組系數n即為輪胎數量。

如3.2 節討論，道面積水分布可改變輪組內部滑水行為，影響輪組阻力特征。對此，將單個輪胎阻力結果疊加得出輪組阻力曲線，并與式（2）計算位移阻力理論解求得比值，得出不同加載荷位時輪組系數n隨滑行速度的變化規律，如圖11 所示。

圖11 輪組系數隨滑行速度的變化曲線Fig.11 Variation cures of wheel configuration coefficient versus taxing speed

如圖11（a）所示，分別以左輪、右輪和輪組平均阻力理論解為底，計算輪組系數曲線并進行比較，可以看出當輪組作用于A2 加載荷位時，積水較厚的左輪阻力理論解大于右輪，因而圖11（a）中左輪曲線位于右輪曲線下方，曲線變形規律基本一致；由平均阻力曲線計算輪組系數中位數等于1.87，略低于期望值2.0。

如圖11（b）所示，在飛機著陸滑行至脫離跑道速度區間內（150～250 km/h），平整道面條件下輪組系數曲線在2.0 上下波動，偏離度不超過2%，輪組系數簡化方式基本適用。對比A1 加載不同荷位結果，輪組系數曲線整體在2.0 以下范圍波動，實際輪組積水阻力低于式（2）計算理論值?？紤]到飛機自動剎車系統，防滑系統和自動擾流系統均以輪胎滑行速度作為激活條件，當輪胎滑行速度未能降低到規定值時，上述3 類措施無法及時發揮減速作用，導致飛機著陸滑跑距離延長，影響飛機著陸安全。

4 結 論

1）輪胎發生滑水前道面積水附加阻力隨滑行速度增加而逐步增大，峰值時刻對應輪胎臨界滑水速度，與NASA 公式結果相差在5%以下；發生滑水后積水阻力非線性減小，持續影響飛機滑行狀態，較道面支撐力更適合作為滑水分析指標。

2）飛機著陸時輪胎與積水道面高速接觸可發生瞬時滑水，后隨減速滑行脫離滑水狀態，影響道面流域特征；相同參數條件下著陸滑行道面積水阻力始終低于起飛滑行過程，著陸階段臨界滑水速度較起飛階段低8.3%～10.6%，表明飛機著陸過程滑水風險更高，符合事故統計規律。

3）輪轍變形可改變道面積水條件，飛機輪組滑水過程存在時空分布差異，空間上表現為輪組前緣動水壓強非對稱分布規律，時間上表現為左右輪胎并非同時達到臨界滑水狀態，vP差值隨積水厚度差增加而增大。

4）僅對平整道面理想積水條件，道面積水阻力輪組系數可近似按輪胎數量計算；有輪轍道面條件下輪組積水阻力低于理論值，輪組系數中位數低于2.0，可拖延減速滑行過程，對預防沖出跑道事故較為不利。