全鋼子午線輪胎溫度場仿真研究

熊春明1，臧孟炎1，晏 寧2

(1.華南理工大學機械與汽車工程學院，廣東 廣州 510640)(2.萬力輪胎股份有限公司，廣東 廣州 511400)

摘要：針對輪胎運動過程中的生熱問題，建立295/80R22.5型全鋼子午線輪胎溫度場有限元模型。使用ABAQUS的穩態傳輸功能，分析輪胎在80km/h速度下的穩態滾動，獲得輪胎模型中各單元的應力應變；使用MATLAB編寫程序計算橡膠單元生熱率，編寫ABAQUS用戶子程序為輪胎溫度場模型定義熱源；為溫度場模型定義熱邊界條件，仿真計算輪胎溫度場。結果表明：在80km/h工況下，295/80R22.5型全鋼子午線輪胎的胎面和胎側溫度與實驗結果基本一致；仿真結果得出高溫區域出現在胎面和胎肩部位，適當增加2號帶束層和3號帶束層間帶束膠的厚度能降低輪胎最高溫度。

關鍵詞：全鋼子午線輪胎；生熱；溫度場；有限元分析

車輛行駛過程中，輪胎受交變載荷作用。由于橡膠材料存在彈性遲滯效應，從而產生大量的熱，導致輪胎溫度升高。溫度升高使得橡膠材料軟化，輪胎剛度變低，變形增大，能量損耗增加；同時，溫度過高，橡膠材料會發生化學老化與過熱破壞，加快輪胎的疲勞破壞，導致輪胎過早發生疲勞失效[1]。因此，全鋼子午線輪胎的溫度場問題受到國內外學者和輪胎企業的重視。

傳統的實驗方法是研究輪胎溫度場的重要手段[2-3]。隨著計算機性能的提高，仿真分析軟件的不斷完善，有限元方法越來越多地用于輪胎溫度場的研究。EBBOTT等[4]將輪胎與外界的熱交換定義為熱對流，使用固定熱對流系數的方法對光面輪胎溫度場進行預測。吳福麒[5]使用非中心對稱的局部線性化模型對膠料的損耗特性進行表征，并將光滑凸旋轉曲面的對流換熱關系用于研究輪胎換熱系數，研究了滾動輪胎溫度場。王學瑞等[6]使用專用有限元軟件TYSYS研究了額定負荷及超負荷條件下的輪胎溫度場。

輪胎溫度場非常復雜，主要表現在以下2個方面： 1)輪胎的溫升問題涉及復雜的熱力耦合過程；2)滾動輪胎表面的熱邊界條件難以確定。本文以輪胎溫度場的單向解耦思想[7]為理論基礎，以295/80R22.5型全鋼子午線輪胎為研究對象，以ABAQUS軟件為工具，建立了輪胎溫度場分析模型。復雜的熱力耦合過程由3個部分構成：穩態滾動分析、生熱率計算和熱傳導分析。熱邊界條件使用反演法[8]進行確定，即根據實驗測試所得的溫度值，不斷調整換熱系數，直到輪胎表面溫度的仿真值和實驗值基本一致。前人所做的工作中，多使用常溫下的損耗因子，沒有考慮到損耗因子會隨溫度的改變而變化。本文根據仿真分析得到的溫度結果，將最初使用的損耗因子更新為仿真模型中對應位置溫度下的損耗因子，對溫度結果進行迭代，最后得到輪胎溫度場。

1 建立溫度場仿真模型

溫度場的建模過程主要分為3個部分：建立3D輪胎有限元模型、計算生熱率和分析熱傳導。具體過程為：首先從有限元模型中輸出輪胎各位置的應力應變，然后以應力應變為基礎計算各單元生熱率，最后將生熱率導入熱傳導模型中，添加熱邊界條件，得到輪胎溫度場模型，如圖1所示。

1.1 輪胎三維有限元模型

全鋼子午線輪胎多用于大型客車，主要由胎面、胎側、胎圈及胎體組成。輪胎主要由具有超彈性特性的橡膠材料以及起支撐作用的骨架材料(橡膠-簾線復合材料)構成。

建模過程中，首先對輪胎二維截面劃分網格，然后使用ABAQUS建立二維有限元模型及三維有限元模型。考慮到輪胎使用過程中存在大變形情況，橡膠材料采用Yeoh模型進行模擬；為了更好地輸出簾線力以及方便劃分網格，橡膠-簾線復合材料采用rebar單元進行模擬[9]。輪輞及路面采用分析剛體進行模擬。為了提高計算效率，輪輞安裝和充氣兩個步驟在二維模型中完成。得到的輪胎三維有限元模型如圖2所示。

進行輪胎徑向剛度有限元仿真和實驗，結果如圖3所示。從圖中可以看出,仿真曲線和實驗結果趨勢一致，斜率也很接近，驗證了有限元模型建立的正確性。

1.2 生熱率計算

輪胎在滾動過程中，輪胎內部存在周期性的應力和應變。橡膠屬于黏彈性材料，變形過程中應力應變不同步，存在相位差δ，形成損耗角，如圖4(a)所示。循環載荷作用下，應力應變的不同步會導致材料應力應變曲線出現一個滯后環，如圖4(b)所示，這個滯后環所包含的面積即為橡膠材料的滯后損失。若橡膠材料承受正弦形態的應力應變，則該應力應變可用數學函數表示為：

ε=ε0sin(ωt)

(1)

σ=σ0sin(ωt+δ)

(2)

式中：ε為應變；σ為應力；ε0為應變幅值；σ0為應力幅值；ω為角頻率；δ為滯后角，t為時間。

根據式(1)和式(2)的關系，一個周期T內單位體積的能量損耗為：

(3)

將式(1)和式(2)代入式(3)可得：

E=πσ0 ε0sinδ

(4)

對于線性黏彈性材料，δ與ε0和σ0無關，可直接按照式(4)計算能量損耗。對于非線性黏彈性材料，理論上不存在與幅值無關的滯后角δ。本文按照文獻[8]中的損耗因子處理方法，采用等效損耗因子替代損耗因子，式中sinδ用等效黏性損耗系數

代替。即：

(5)

式(5)中的

不僅與ω相關，還與σ0,ε0相關。

根據彈性理論，材料的變形可用等效應力應變表示。等效應變為：

(6)

等效應力為：

(7)

以上能量損耗計算方法只適用于應力應變為正弦曲線的情況，然而輪胎在滾動過程中，其應力應變并非以正弦曲線的形式變化。對于非正弦變化的應力應變，通過傅里葉變換[10]可得到計算能量損耗所需的應力幅值σ0和應變幅值ε0。變換公式如下：

(8)

(9)

應變幅值為：

(10)

應力幅值為：

(11)

式(8)、(9)中m取30即可很好地擬合應力應變。

根據式(4)、(10)、(11)計算單位體積橡膠一個周期內的能量損耗：

(12)

按照上述方法，使用MATLAB編寫程序計算橡膠單元的能量損耗：首先從穩態滾動模型中提取單元應力應變，再將應力應變轉化為等效應力應變；然后對等效應力應變作傅里葉變換得到等效應力應變幅值；最后根據式(12)計算能量損耗。損耗因子隨溫度的升高有所降低，初次計算能量損耗時選擇室溫25℃時各種膠料的損耗因子。

1.3 溫度場建模

自由滾動狀態下，輪胎滾動一周，輪胎周向方向每個截面所經歷的應力應變狀態是相同的，可以認為輪胎周向不存在溫度差。實驗表明，輪胎運行2h后生熱與散熱達到平衡，溫度場達到穩態狀態。故輪胎斷面溫度場可代替整個輪胎的溫度分布。基于上述分析，本文將輪胎溫度場模型簡化為二維模型。溫度場模型的建立主要分為以下2個部分：定義熱源、定義熱邊界條件。

1.3.1 定義熱源

ABAQUS軟件中不能直接在CAE模塊或者inp文件中定義熱源，而是需要使用子程序進行定義。將計算得到的生熱率數據，按單元編號順序排列，以備定義熱源的子程序調用。本文在定義內熱源的過程中用到3個子程序：子程序UEXTERNALDB用于將生熱率數據讀取到程序中以備使用，子程序HETVAL直接用于定義熱源，但其不能識別單元編號，需要調用子程序SDVINI進行單元編號識別[11]。

1.3.2 熱邊界條件

輪胎與外界的熱邊界包括與地面及空氣接觸的胎面邊界、與空氣接觸的胎側邊界、與輪輞接觸的胎圈邊界以及與空氣接觸的內腔邊界。在輪胎與外界的熱交換中，熱對流為輪胎與外界熱交換的主要形式[1]。實際運動過程中，胎圈和輪輞間主要以熱傳導形式進行熱交換；胎面與路面間主要以熱傳導形式進行熱交換，其余表面直接和空氣進行熱交換，屬于熱對流。為簡化溫度場模型，將胎圈與輪輞間的熱交換、胎面與路面間的熱交換都簡化為熱對流，通過調整對流換熱系數的方式定義不同位置的熱邊界條件。

本文采用反演法[8]確定圖5所示的胎面、胎側、胎圈及內腔的熱對流系數。該方法根據實驗實測溫度不斷調節對流換熱系數，最終使溫度場模型表面溫度和實測輪胎表面溫度分布達到一致。

測量輪胎溫度時，實驗室溫度為25℃，因此設置與胎面、胎側及胎圈接觸的空氣溫度為25℃，輪胎初始溫度為25℃。為了確定輪胎內表面與內腔空氣的熱邊界條件，還需要確定內腔空氣溫度。輪胎運行過程中內部生熱將導致內腔溫度升高，在輪胎溫度達到平衡狀態后，內腔溫度與氣壓也會達到平衡狀態。由于輪胎內腔體積變化小，可將這一過程視為理想氣體等容變化。由理想氣體狀態方程(式(13))可以根據不同狀態下的氣壓估算輪胎內腔在不同狀態下的溫度。

(13)

式中：P1為輪胎運行前內腔氣壓，實測值為0.85MPa；P2為輪胎運行至動態平衡后內腔氣壓，實測值為1.39MPa；T1為輪胎運行前內腔溫度，實測值為25℃；T2為輪胎運行至動態平衡后內腔溫度，根據式(13)，代入數據，計算結果為41℃。

2 溫度場仿真計算與評價

2.1 溫度場仿真計算

按上述邊界條件計算所得到的295/80R22.5輪胎模型的溫度場如圖6所示，胎肩及胎面位置溫度較高，最高溫度出現在胎肩的帶束膠位置，為135.4℃。

由于損耗因子隨溫度的升高而有所降低，而前文計算生熱率時使用了25℃的損耗因子，不能正確反映輪胎各部位的實際能量損耗，因此需要根據輪胎不同部位的溫度重新選擇損耗因子，計算生熱率。具體操作方法如下：首先從295/80R22.5模型結果文件中輸出不同部位每個單元的溫度數據，接著計算每個部位單元溫度的平均值，最后根據各個部位單元溫度的平均值選擇對應溫度下的損耗因子重新計算生熱率。更新損耗因子后的溫度場如圖7所示，胎肩及胎面位置仍然為高溫區，最高溫度出現在胎肩帶束膠位置，降為125.6℃。將這一過程反復實施，可以不斷提高輪胎溫度場的仿真精度，直到兩次計算結果的最高溫度差值小于0.2℃為止。

295/80R22.5修改帶束模型溫度場計算結果如圖8所示，最高溫度為122.4℃，相比295/80R22.5原始模型有所降低，說明修改帶束結構可降低胎肩溫度。

2.2 溫度場仿真結果評價

為了對溫度場模型進行驗證，使用手持式紅外測溫儀對輪胎外表面進行溫度測試。溫度測試位置如圖9所示，仿真與實驗的溫度對比結果見表1。

由表1可知，大多數位置仿真和實驗的誤差都在5%以內，但tread1位置誤差達11.4%。由于輪胎結構左右對稱，載荷左右對稱，故輪胎胎面溫度也應該大致左右對稱。然而，tread1位置的實驗測試溫度為42.4℃，tread5位置的實驗測試溫度為47.3℃，并且42.4℃和胎面其他位置溫度存在明顯差異，因此可以認為tread1位置溫度的實驗值明顯低于仿真值是實驗誤差所致。溫度場仿真分析結果與實驗的一致性，說明了仿真分析方法的有效性。

因為最高溫度出現在胎肩部位的帶束膠位置，所以對帶束端部結構的改進有望降低輪胎斷面最高溫度。圖8所示為適當增加2號帶束與3號帶束間帶束膠的厚度后輪胎斷面的溫度分布，盡管最高溫度發生的位置不變但是最高溫度下降了3.2℃，說明了帶束端部的結構優化對降低溫度以及提高輪胎疲勞壽命的重要性。

3 結束語

本文從輪胎穩態滾動有限元分析開始，通過穩態滾動仿真、生熱率計算和輪胎熱傳導分析得到輪胎溫度場仿真分析有限元模型。在此基礎上，根據輪胎仿真溫度結果對橡膠的損耗因子進行更新，重新計算溫度場，得到的結果相比前人的方法更接近實際情況。仿真結果與實驗的一致性也說明了仿真分析方法的有效性。本文開發的輪胎溫度場仿真模型可以為輪胎設計人員提供針對性的建議，在輪胎的開發過程中對易發生破壞的位置采取對應措施，從而縮短產品開發時間。下一步，筆者將以此為基礎，深入進行輪胎疲勞壽命分析和優化研究工作。

參考文獻：

[1] 張士齊. 輪胎力學與熱學[M]．北京：化學工業出版社，1988:280-283.

[2] 黃夢溪, 李子然, 夏源明. 滾動輪胎表面紅外成像測溫的偏差及其修正[J]．實驗力學, 2013, 28(5):595-600.

[3] 李勇, 左曙光, 段向雷,等. 基于試驗的輪胎溫度場分布及影響因素分析[J]．同濟大學學報(自然科學版), 2012, 40(8):1249-1253.

[4] EBBOTT T G, HOHMAN R L, JEUSETTE J P, et al. Tire temperature and rolling resistance prediction with finite element analysis[J]．Tire Science and Technology, 1999, 27(1): 2-21.

[5] 吳福麒. 輪胎穩態滾動溫度場的有限元分析[D]．合肥：中國科學技術大學, 2009.

[6] 王學瑞, 王澤君, 王友善,等. 無內胎全鋼輕型載重子午線輪胎穩態滾動溫度場有限元分析[J]．橡膠工業, 2012, 59(3):173-176.

[7] SAMUEL K C, RICHARD N D. Heat generation in aircraft tires under free rolling conditions[M].Washington：NASA Langley Technical Report Server, 1982.

[8] 王曉軍. 子午線輪胎溫度場有限元分析與測試[D]．合肥：中國科學技術大學, 2004.

[9] 熊春明, 臧孟炎, 周濤,等. 全鋼子午線輪胎疲勞破壞現象仿真與優化設計研究[J]．輪胎工業, 2016, 36(9):515-519.

[10] SHIDA Z, KOISHI M, KOGURE T, et al. A rolling resistance simulation of tires using static finite element analysis[J]．Tire Science and Technology, 1999, 27(2): 84-105.

[11] 陳秋紅. 載重子午線輪胎穩態溫度場分析[D]．上海：東華大學, 2008.

Temperature field simulation research on all-steel radial tire

XIONG Chunming1,ZANG Mengyan1,YAN Ning2

(1.School of Mechanical and Automotive Engineering, South China University of Technology,Guangdong Guangzhou 510640，China) (2.Wanli Tire Corporation Limited,Guangdong Guangzhou 511400,China)

Abstract：Aimed at tire heating problem during rotation, temperature field model of all-steel radial tire 295/80R22.5 is built. Based on ABAQUS, steady-state rolling model of tire under a speed of 80km/h is built to obtain strains and stresses of each element. With the help of MATLAB, heating rate of rubber elements are calculated according to related equations. ABAQUS subroutines are programmed to define heat source of tire temperature field model. And then, heat boundary conditions are added to temperature field model. Tire temperature field are finished. The results indicate that, under a speed of 80km/h, simulation temperature distribution of tread and side wall is the same as experimental results. Temperatures of tread and tire shoulder are higher than other positions. What's more, appropriately increasing thickness of the composite between 2# belt and 3# belt can lower the temperature effectively.

Key words：all-steel radial tire; tire heating; temperature field model; finite element analysis

DOI：10.3969/j.issn.2095-509X.2017.02.007

收稿日期：2016-12-28

作者簡介：熊春明(1992—)，男，重慶人，華南理工大學碩士研究生，主要從事汽車CAE方向的研究。

中圖分類號：U463.341

文獻標識碼：A

文章編號：2095-509X(2017)02-0035-05