液力緩速器部分充液流場大渦模擬及特性預測

李雪松，程秀生，苗麗穎

(吉林大學 汽車工程學院，吉林 長春，130022)

液力緩速器作為一種輔助制動裝置，依靠工作輪內液流的作用將車輛的動能轉化為液體的熱能，再通過冷卻器散熱的方式實現車輛制動[1]。當液力緩速器工作在部分充液工況下，其內部是復雜的黏性、不可壓縮、非穩態氣-液兩相流動，流場測試試驗較復雜，且設備昂貴。隨著計算機技術和計算流體動力學的迅速發展，對液力緩速器內流場進行數值模擬的方法比較經濟可行。為全面揭示液力緩速器流場的內部特性以及內外特性之間的變化規律，目前國內學者常用雷諾平均法(RNS)進行流場的數值計算[2-4]，其弊端是忽略了湍流的若干細節，所以，對封閉的雷諾時均方程進行求解并不能有效地模擬流動中分離、渦旋以及擴散等現象[5]。本文作者借助CFD軟件平臺，采用更為準確的大渦模擬法對液力緩速器部分充液流場進行數值模擬，獲取液力緩速器內部流場速度和壓力分布特性，以及不同充液率下內流場結構的變化和兩相體積分布情況，并基于流場數值解計算液力緩速器的外特性。

1 控制方程及計算模型

大渦模擬(LES)是介于直接數值模擬(DNS)與雷諾平均法(RANS)之間的一種湍流數值模擬方法,是目前CFD研究的熱點之一[6-8]，其基本假設是：動量、能量、質量方程及其他標量主要由大渦輸運；流動的幾何和邊界條件決定了大渦的特性，而流動特性主要在大渦中體現；小尺度渦旋受幾何合邊界條件影響較小，并且各向同性，直接求解大渦，而小尺度渦對大尺度渦運動的影響則通過一定的模型在針對大尺度渦的瞬時N-S方程體現出來[9]。

LES的控制方程是對N-S方程在波數空間或者物理空間進行過濾得到的，通過去掉比過濾寬度或者給定物理寬度小的渦旋，得到其控制方程：

為使控制方程組封閉，必須建立亞格子尺度模型(SGS模型)，目前應用最廣泛的漩渦黏性模型方程為：

式中：δij是“Kronecker delta”符號(當i=j時，δij=1；當 i≠j時為亞格子尺度的湍動黏度，根據Smagorisnsky-Lily模型定義，其中LS為網格濾波寬度。

式中：k為 VonKarman常數，k＝0.418 7；CS為Samagorin常數，CS=0.1；d為到最近的壁面的距離；V為計算單元體積。

根據液力緩速器樣機建立的三維模型以及提取的計算區域網格模型分別如圖1和圖2所示。定子和轉子循環圓有效直徑分別為293 mm和296 mm，其葉片數量分別為34和36，均為40°前傾葉片。考慮到液力緩速器定子內有進出油口，定子葉片厚度不同，工作流道幾何結構不是周期對稱，因此，選取整個流道空間作為計算模型，網格單元總數為794 640個。

圖1 液力緩速器CAD模型Fig.1 Hydraulic retarder CAD model

圖2 液力緩速器計算網格模型Fig.2 Computational mesh models of hydrodynamic retarder

2 數值計算方法

2.1 控制方程的分離式求解法

流場數值計算的本質就是求解離散后的控制方程組，其基本過程是在空間上將計算域離散成許多小的體積單元，在每個體積單元上對離散后的控制方程組進行求解。求解方法可分為分離式求解法和耦合式求解法。分離式求解法通常用于不可壓縮流體數值計算，其求解時同時考慮所有控制單元，順序、逐個地求解各變量代數方程組。目前工程上使用最廣泛的壓力修正法，其單個時間步的求解過程如圖3所示。本文運用大渦模擬法對液力緩速器內不可壓縮的三維瞬態流場進行模擬，也必須采用分離式求解法求解其控制方程組[10-11]。

圖3 分離求解法計算流程圖Fig.3 Flow chart for segregated method

2.2 可動區域耦合計算的滑動網格法

液力緩速器的瞬態特性主要體現在其工作腔內液體的流動狀態為湍流非穩態流動，當轉子轉速一定、工作腔內充液率保持不變時，轉子和定子之間相對位置的不斷變化引起工作流道內的流場特性呈現周期性的變化，當緩速器工作腔內的充液率及轉子轉速不斷變化時，工作流道內呈現復雜的非穩態流動形態。

為了對轉子和定子工作流道進行統一計算，對轉子和定子的接合面采用了滑動網格理論來傳遞兩個區域之間的參數。滑動網格法是模擬多移動參考系流場的最精確的瞬態計算方法，其基本原理為：將計算區域中相鄰子單元的分界面互相耦合形成“網格分界面”，數值計算中單元區域在離散步驟中沿著網格分界面相互滑動(旋轉或平移)，同時分界面也隨時間發生變化(如圖 4所示)[12]。滑動網格模擬的瞬態問題大部分是時間周期性的，即計算區域的速度是周期復現的。

設T為瞬態計算的周期，在計算區域的一些流動特性函數Φ為

轉子在一個周期T內運動2π，又其轉速為1 200 r/min，則計算周期T可用下式計算：

圖4 網格滑移前后葉輪之間相對位置Fig.4 Relative position of impellers for sliding mesh theory

單位時間步長以不大于單位網格單元大小為選取原則，文中選取的單位時間步t=0.12 ms。

液力緩速器工作在部分充液狀態時，氣相與液相之間相互混摻，將混合模型與歐拉模型交替運用在其兩相流數值計算中，在多相流模型中首先采用混合模型獲得初始解，并以此作為起點，啟動歐拉模型繼續計算。在制動特性研究中往往關心液力緩速器內外特性的關系，即轉子轉速、充液率以及制動轉矩的變化規律，而這主要是由液相起支配作用，因此設定液相為主相，氣相為附加相。

假設液力緩速器內的循環液體作等溫流動且沒有泄露，認為其在流動過程中為不可壓縮黏性流體，且流體與葉輪間的流固耦合作用不致引起流道的變形。流道內壁與葉片表面近壁處速度場計算采用速度無滑移邊界條件，進出油口采用壓力邊界條件。

3 數值計算結果分析

為了準確掌握液力緩速器內部分充液時的流場分布特性，對轉子轉速為1 200 r/min，充液率f分別為25%和75%的數值計算結果進行對比分析。

3.1 節面流場特性分析

液力緩速器內不同充液率下葉片壓力面的靜壓力分布如圖5所示。較高的充液率能導致較大的靜壓力，因為當轉子轉速相同時，充液率越大，流體密度也越大，則引起的靜壓力越大。將不同流動單元中壓力面的靜壓分布(圖 5中(Ⅰ)～(Ⅳ))對比分析發現：其分布趨勢相同，即從循環圓中心到外環，靜壓力遞增，入口處出現局部高壓區，主要是由于定-轉子之間過大的轉速差導致此處所受液流沖擊最大所致；在循環圓中心，渦旋和二次流現象頻繁，因此存在環形低壓區。充液率較高時，葉片內的壓差較大，液力損失也大。

不同充液率下吸力面的靜壓力分布如圖6所示。由于入油口分布在不同流動單元內葉片吸力面上，因此不同流動單元內吸力面靜壓分布差別較大(圖 6中(Ⅰ)～(Ⅳ))。轉子葉片中，液流離心力作用使得葉片根部出現高壓區；定子吸力面內沒有大面積高壓區，低壓區則是由于沖向壓力面的液流在吸力面產生了尾流造成的。由圖6中(Ⅲ)和(Ⅳ)可見：定子葉片的進油口附近壓力較低，正是為了保證液流以高速進入工作腔。

3.2 弦面流場特性分析

圖7所示為不同弦面相對速度矢量分布。中間弦面的渦旋及二次流動現象幾乎均布在各流動單元內兩葉輪的交接處，這主要是由于壓力面和吸力面之間存在較大壓差。此外，該弦面內流體質點間的相對劇烈運動也是導致渦旋形成的重要因素，充液率越高，渦旋在各流動單元內分布越均勻[13-14]。如圖中A處所示，充液率較低時，在進油口較集中的部分流動單元內出現了逆流以及小渦旋。

圖8所示為相對應的弦面靜壓力分布圖。從整體來看，不同充液率時壓力分布趨勢相同。充液率為75%時，整個弦面內的壓力差以及吸力面與壓力面之間的壓差較25%充液時大。壓力差是使液流偏離循環流動方向形成二次流的重要因素[15]，與圖7對比可發現逆流及二次流較明顯的部分流動單元內，其靜壓值也偏低。

圖5 壓力面靜壓分布Fig.5 Static pressure distribution of pressure surface

圖6 吸力面靜壓分布Fig.6 Static pressure distribution of suction surface

圖7 弦面相對速度矢量分布Fig.7 Velocity vector distribution of chord surface

圖8 弦面靜壓分布Fig.8 Static pressure distribution of chord surface

3.3 液相體積分布

圖9所示為不同充液率時的液相體積分布。由于液體在工作腔內的循環流動，使液相主要集中在流道外環，而氣相分布在內環。充液率為75%時，工作流道的大部分區域被液相占據，而在循環圓中心區域形成氣環；充液率為25%時，氣相占據了工作區域的大部分空間。

4 性能預測與試驗對比

通過對液力緩速器內流場進行瞬態數值計算，可以得到轉子葉片表面的速度和壓力梯度，據此獲得轉子葉片表面各單元的壓力和剪切力，再通過對旋轉軸求矩可得到該單元對旋轉軸的作用力矩。對所有轉子葉片表面單元進行積分求和即可得到流體作用在轉子上的總力矩。滑動網格模擬的瞬態問題與時間t有關，經過若干周期，轉矩T已經隨時間進行周期性變化，因此可以進行數據采集。

表1和表2所示分別為數值計算結果以及在汽車動態模擬國家實驗室液力傳動試驗臺獲得的試驗數據。由表1，2和圖10可知：基于三維瞬態流場數值解的轉矩計算結果與試驗測試結果吻合較好，最大相對誤差在8%以內，說明采用的數值計算方法可以較準確地預測液力緩速器的性能。

圖9 液相體積分布Fig.9 Volume distribution of liquid-phase

表1 數值計算結果Table 1 Numerical simulation results

表2 試驗數據結果Table 2 Experimental data results

圖10 數值計算結果與試驗結果對比Fig.10 Comparison of numerical results with experiments

5 結論

(1) 通過定性對比發現，采用大渦模擬法，并利用滑動網格理論處理定-轉子之間的流動參數傳遞，可以較準確地模擬液力緩速器內氣-液兩相流場的真實流動狀態，并預測其性能。

(2) 外特性數值計算結果與試驗結果吻合良好，因此可以用軟件計算來代替部分試驗工作，其結果可用于液力緩速器的設計及其結構優化。

[1]何仁. 汽車輔助制動裝置[M]. 北京: 化學工業出版社, 2005:19-21.HE Ren. Automotive auxiliary braking system[M]. Beijing:Chemical Industry Press, 2005: 19-21.

[2]何延東, 馬文星, 劉春寶. 液力偶合器部分充液流場數值模擬與特性計算[J]. 農業機械學報, 2009, 40(5): 24-28.HE Yan-dong, MA Wen-xing, LIU Chun-bao. Numerical simulation and characteristic calculation of hydrodynamic coupling[J]. Transaction of CSAM, 2009, 40(5): 24-28.

[3]嚴軍, 何仁, 魯明. 液力緩速器變葉片數的三維數值模擬[J].江蘇大學學報: 自然科學版, 2009, 30(1): 27-31.YAN Jun, HE Ren, LU Ming. Numerical simulation of hydraulic retarder with different blade number[J]. Journal of Jiangsu University: Natural Science Edition, 2009, 30(1): 27-31.

[4]張德生, 趙繼云, 劉立寶, 等. 基于 CFD 的桃形腔偶合器流場分析及結構優化[J]. 中國礦業大學學報, 2010, 39(5):687-692.ZHANG De-sheng, ZHAO Ji-yun, LIU Li-bao, et al. Flow field analysis and structure optimization of peach shaped chamber hydrodynamic coupling based on CFD[J]. Journal o f China University of Mining & Technology, 2010, 39(5): 687-692.

[5]才委, 馬文星, 褚亞旭, 等. 液力變矩器導輪流場數值模擬與試驗[J]. 農業機械學報, 2007, 38(8): 11-14.CAI Wei, MA Wen-xing, CHU Ya-xu, et al. Numerical simulation and experimental research on flow field in the stator of a torque converter[J]. Transactions of the Chinese Society for Agricultural Machinery, 2007, 38(8): 11-14.

[6]Brandt T T. Usability of explicit filtering in large eddy simulation with a low-order numerical scheme and different subgrid-scale models[J]. International Journal for Numerical Methods in Fluids, 2008, 57(7): 905-928.

[7]ZHANG Yan-hong, YANG Chao, MAO Zai-sha. Large eddy simulation of the gas-liquid flow in a stirred tank[J]. AICHE Journal, 2008, 54(8): 1963-1974.

[8]Grigoriadis D G E, Bartzis J B, Goulas A. Efficient treatment of complex geometries for large eddy simulations of turbulent flows[J]. Computers and Fluid, 2004, 33(2): 201-222.

[9]王福軍. 計算流體動力學分析—CFD軟件原理與應用[M]. 北京: 清華大學出版社, 2004: 139-142.WANG Fu-jun. Computational fluid dynamics analysis:Principle and application of CFD software[M]. Beijing: Tsinghua University Press, 2004: 139-142.

[10]Felten F, Fautrelle Y, Du Terrail Y, et al. Numerical modeling of electrognetically-riven turbulent flows using LES methods[J].Applied Mathematical Modelling, 2004, 28(1): 15-27.

[11]Feiz A A, Ould-Rouis M, Lauriat G. Large eddy simulation of turbulent flow in a rotating pipe[J]. International Journal of Heat and Fluid Flow, 2003, 243(3): 412-420.

[12]FLUENT Inc. FLUENT user’s guide[M]. Fluent Inc, 2003:239-244.

[13]Mary I, Sagaut P. Large eddy simulation of flow around an airfoil near stall[J]. AIAA Journal, 2002, 40(6): 1139-1145.

[14]Julian R E, Smolarkiewicz K. Eddy resolving simulations of turbulent solar convection[J]. International Journal for Numerical Method in Fluids, 2002, 39(9): 855-864.

[15]Bai L, Fieblg M, Mitra N K. Numerical analysis of turbulent flow in fluid couplings[J]. ASME Journal of Fluids Engineering,1997, 119(3): 569-576.