基于STAR-CCM+的某發(fā)動(dòng)機(jī)冷卻風(fēng)扇數(shù)值模擬

唐 琳，譚禮斌，黃 燦，劉小強(qiáng)，余千英，何 丹，冷小麗

（隆鑫通用動(dòng)力股份有限公司技術(shù)中心基礎(chǔ)研究所CFD研究室，重慶400039）

0 引言

發(fā)動(dòng)機(jī)冷卻風(fēng)扇主要作用是利用空氣的流動(dòng)，對(duì)運(yùn)行中的發(fā)動(dòng)機(jī)進(jìn)行冷卻[1]。冷卻風(fēng)扇的風(fēng)量是評(píng)估風(fēng)扇性能的一個(gè)重要指標(biāo)。隨著CFD（Computational Fluid Dynamics）仿真分析技術(shù)在工業(yè)領(lǐng)域的廣泛運(yùn)用和發(fā)展，采用數(shù)值仿真方法對(duì)旋轉(zhuǎn)機(jī)械內(nèi)部流場(chǎng)進(jìn)行可視化分析，已逐漸成為整個(gè)行業(yè)的趨勢(shì)[2]。本文基于計(jì)算流體動(dòng)力學(xué)基本理論，采用通用流體分析軟件STAR-CCM+對(duì)某發(fā)動(dòng)機(jī)冷卻風(fēng)扇進(jìn)行穩(wěn)態(tài)流場(chǎng)特性分析，預(yù)測(cè)其在給定實(shí)驗(yàn)測(cè)試工況條件下的風(fēng)量大小，并與實(shí)驗(yàn)測(cè)試的風(fēng)扇風(fēng)量值對(duì)比，驗(yàn)證CFD仿真分析的可靠性，為后續(xù)冷卻風(fēng)扇的設(shè)計(jì)及結(jié)構(gòu)改進(jìn)提供相應(yīng)的理論基礎(chǔ)和分析支撐。

1 物理模型

某發(fā)動(dòng)機(jī)冷卻風(fēng)扇的實(shí)物圖如圖1a所示，采用Unigraphics NX 8.0軟件按照冷卻風(fēng)扇實(shí)物以1∶1繪制得出3D模型（圖1b）。將三維模型導(dǎo)入流體分析軟件STAR-CCM+11.06后，依據(jù)實(shí)驗(yàn)測(cè)試模型建立計(jì)算域模型。采用多面體網(wǎng)格和邊界層網(wǎng)格類型對(duì)風(fēng)扇計(jì)算域模型進(jìn)行網(wǎng)格劃分，獲得相應(yīng)的計(jì)算域網(wǎng)格模型圖（圖1c）。本文對(duì)網(wǎng)格數(shù)量和球面入口尺寸進(jìn)行網(wǎng)格無關(guān)性研究，對(duì)過渡區(qū)域、旋轉(zhuǎn)域、非旋轉(zhuǎn)域進(jìn)行體加密、葉尖加密，研究網(wǎng)格數(shù)量對(duì)風(fēng)量計(jì)算結(jié)果的影響；改變球面尺寸，探究其對(duì)風(fēng)扇風(fēng)量計(jì)算結(jié)果的影響。體加密示意圖及球面入口示意如圖2所示。

圖1 冷卻風(fēng)扇模型示意圖

圖2 網(wǎng)格加密示意圖及球面入口示意圖

2 數(shù)學(xué)模型

流體流動(dòng)滿足質(zhì)量守恒定律、動(dòng)量守恒定律和能量守恒定律。本文數(shù)值求解中不涉及溫度場(chǎng)，因此流動(dòng)控制方程如下：

（1）連續(xù)方程[3]

該方程表述為在同一時(shí)間內(nèi)某一封閉流體控制體表面的體積流量為零，即流入的體積流量等于流出的體積流量。不可壓縮流體下，密度ρ為常數(shù)的質(zhì)量守恒方程如式（1）所示。

（2）動(dòng)量方程（N-S方程）[4]

該方程可表述為某一流動(dòng)控制體，其動(dòng)量的積累速率等于作用于控制體上的外力的總和。

式中：ui、uj是平均速度分量，xi、xj為坐標(biāo)分量，p是流體微元體上的壓力，μeff是湍流有效黏性系數(shù)。

本文選用STAR-CCM+中提供的Realizable k-ε湍流模型進(jìn)行數(shù)值計(jì)算，因此在數(shù)值求解計(jì)算中需要求解k-ε湍流模型方程[4]：

式中：Gk為速度梯度產(chǎn)生的湍動(dòng)能項(xiàng)，Gb為浮力產(chǎn)生的湍動(dòng)能項(xiàng)，YM表為脈動(dòng)擴(kuò)張項(xiàng)，C1ε，C2ε，C3ε為經(jīng)驗(yàn)常數(shù)，σk，σε分別為與湍動(dòng)能k和耗散率相對(duì)應(yīng)的Prandtl數(shù)，Sk和Sε為用戶自定義的源項(xiàng)。

3 邊界條件

本文風(fēng)扇計(jì)算域模型進(jìn)口采用滯止入口（stagnation inlet）邊界條件，出口采用壓力出口（Pressure outlet）邊界條件，給定的壓力條件為進(jìn)口總壓和出口靜壓，出口壓力值為實(shí)驗(yàn)測(cè)試環(huán)境下各工況的風(fēng)壓值，風(fēng)扇進(jìn)氣溫度25℃.本文采用MRF（Moving Reference Frame）模型實(shí)現(xiàn)動(dòng)—靜面間的數(shù)據(jù)傳遞，給旋轉(zhuǎn)葉片區(qū)域（旋轉(zhuǎn)域）設(shè)定旋轉(zhuǎn)速度和旋轉(zhuǎn)方向，即可實(shí)現(xiàn)風(fēng)扇的運(yùn)轉(zhuǎn)。本文采用穩(wěn)態(tài)求解。

4 模擬結(jié)果分析

4.1 網(wǎng)格無關(guān)性驗(yàn)證

計(jì)算域模型的網(wǎng)格模型劃分方法和網(wǎng)格質(zhì)量對(duì)數(shù)值求解的計(jì)算精度及模擬結(jié)果都具有非常顯著的影響，采用合理的網(wǎng)格參數(shù)控制策略對(duì)關(guān)鍵流動(dòng)區(qū)域進(jìn)行網(wǎng)格細(xì)化控制，對(duì)提高網(wǎng)格質(zhì)量、提高計(jì)算精度都具有非常重要的作用[5]。本文選取19.9 Pa、39.8 Pa兩個(gè)測(cè)試風(fēng)壓進(jìn)行網(wǎng)格無關(guān)性分析，研究不同網(wǎng)格數(shù)量下冷卻風(fēng)扇風(fēng)量的變化情況。圖3為不同網(wǎng)格數(shù)量對(duì)風(fēng)量計(jì)算結(jié)果的影響。圖3表示不同旋轉(zhuǎn)區(qū)域網(wǎng)格尺寸下風(fēng)量計(jì)算結(jié)果的變化。從圖中可以看出，網(wǎng)格數(shù)量達(dá)到600萬(wàn)后，風(fēng)量幾乎沒有變化。同理，圖4和圖5分別表示過渡區(qū)域網(wǎng)格尺寸和非旋轉(zhuǎn)域網(wǎng)格尺寸下風(fēng)量計(jì)算結(jié)果的變化?？傮w上分析，當(dāng)網(wǎng)格數(shù)量達(dá)到600萬(wàn)以上，冷卻風(fēng)扇風(fēng)量的計(jì)算結(jié)果幾乎沒有變化變化。綜上，考慮到實(shí)際計(jì)算能力和葉片變形程度，冷卻風(fēng)扇計(jì)算域模型的網(wǎng)格控制策略如下：非旋轉(zhuǎn)域網(wǎng)格尺寸16 mm，過渡區(qū)域網(wǎng)格尺寸8 mm，旋轉(zhuǎn)域網(wǎng)格尺寸1 mm，葉尖加密尺寸0.25 mm，邊界層層數(shù)6層，邊界層增長(zhǎng)比1.5，風(fēng)扇邊界層厚度0.4 mm，其余邊界層厚度1 mm。網(wǎng)格劃分完成后風(fēng)扇計(jì)算域體網(wǎng)格數(shù)量約為700萬(wàn)。

圖3 不同旋轉(zhuǎn)域網(wǎng)格尺寸對(duì)計(jì)算結(jié)果的影響（非旋轉(zhuǎn)域16 mm；過渡區(qū)域8 mm）

圖4 過渡區(qū)域網(wǎng)格尺寸對(duì)計(jì)算結(jié)果的影響

圖5 非旋轉(zhuǎn)區(qū)域網(wǎng)格尺寸對(duì)計(jì)算結(jié)果的影響

4.2 邊界設(shè)置對(duì)計(jì)算結(jié)果的影響

圖6 為不同氣體入口面對(duì)風(fēng)量計(jì)算結(jié)果的影響。圖6（a）為不同氣體入口面形狀及大小示意圖，圖中2.2倍球面表示球面半徑為風(fēng)扇外輪廓半徑的2.2倍。圖6（b）表示39.8 Pa測(cè)試風(fēng)壓下不同氣體入口面形狀下風(fēng)扇風(fēng)量的計(jì)算結(jié)果，從圖中可以看出，風(fēng)扇罩形狀和立方體形狀的風(fēng)扇風(fēng)量計(jì)算結(jié)果偏大，球面形狀的風(fēng)扇風(fēng)量隨著球面半徑的增大，風(fēng)量逐漸減小，最后趨近于一個(gè)定值。從圖中可以看出，當(dāng)球面半徑大于3.2倍后，風(fēng)扇風(fēng)量值的大小變化較小。因此，本文后續(xù)研究中采用4.2倍球面大小的氣體入口面進(jìn)行冷卻風(fēng)扇風(fēng)量的計(jì)算。

圖6 氣體入口面對(duì)風(fēng)量計(jì)算結(jié)果的影響

STAR-CCM+中的湍流模型有RKE（Realizable KEpsilon）、SKE （Standard K-Epsilon）、Sk-O（Standard K-Omega）、SSTK-O（SST K-Omega）等模型[6]。本文選用4種常用的雙方程模型，分別為SKE模型、RKE模型、SK-O模型和SSTK-O模型，分析湍流模型的選擇對(duì)冷卻風(fēng)扇風(fēng)量計(jì)算結(jié)果的影響。由于不同的湍流模型對(duì)近壁區(qū)網(wǎng)格數(shù)量要求不同，可以用Y+值檢驗(yàn)與壁面最近節(jié)點(diǎn)的位置，確保近壁區(qū)有足夠的節(jié)點(diǎn)數(shù)。在模擬過程中，邊界層網(wǎng)格的Y+值范圍必須滿足相應(yīng)的湍流模型要求，確保數(shù)值模擬的求解精度[6]。如圖7所示為39.8 Pa測(cè)試風(fēng)壓下不同湍流模型對(duì)應(yīng)的邊界層網(wǎng)格Y+值分布圖和冷卻風(fēng)扇風(fēng)量計(jì)算結(jié)果圖。從Y+值分布可以看出，大部分區(qū)域滿足湍流模型需求的Y+值范圍，表明各湍流模型下的邊界層網(wǎng)格Y+值滿足要求。39.8 Pa測(cè)試風(fēng)壓下冷卻風(fēng)扇的風(fēng)量實(shí)驗(yàn)測(cè)試結(jié)果為146.5 m3/h，從風(fēng)量計(jì)算結(jié)果來看，圖中 Realizable K-Epsilon Two-Layer （all Y+）湍流模型下的風(fēng)量計(jì)算結(jié)果為149.7m3/h，與風(fēng)量實(shí)測(cè)值間的誤差最小，表明RKE（ALL Y+）湍流模型對(duì)此冷卻風(fēng)扇的風(fēng)量計(jì)算最適合。

圖7 湍流模型對(duì)風(fēng)量計(jì)算結(jié)果的影響

如圖8所示為39.8 Pa測(cè)試風(fēng)壓下其它因素對(duì)冷卻風(fēng)扇風(fēng)量計(jì)算結(jié)果的影響。從圖中可以看出，采用理想氣體進(jìn)行冷卻風(fēng)扇的風(fēng)量仿真模擬，獲得的風(fēng)量結(jié)果與常密度狀態(tài)的計(jì)算結(jié)果相差不大；采用瞬態(tài)求解器獲得風(fēng)扇風(fēng)量值為148.6 m3/h，與風(fēng)量實(shí)驗(yàn)測(cè)試結(jié)果146.5 m3/h的誤差較小，但CPU計(jì)算時(shí)間消耗約為穩(wěn)態(tài)計(jì)算時(shí)間的4～6倍，計(jì)算成本較高。

圖8 其它因素對(duì)風(fēng)量計(jì)算結(jié)果的影響

4.3 計(jì)算結(jié)果及實(shí)驗(yàn)對(duì)比

采用網(wǎng)格無關(guān)性研究結(jié)論對(duì)冷卻風(fēng)扇計(jì)算域模型進(jìn)行網(wǎng)格劃分，采用實(shí)測(cè)邊界進(jìn)行設(shè)置，湍流模型選擇 Realizable K-Epsilon Two-Layer （all Y+），通過STAR-CCM+11.06流體仿真軟件進(jìn)行穩(wěn)態(tài)數(shù)值求解分析，獲得如圖9所示的冷卻風(fēng)扇Wall Y+圖、速度云圖及速度矢量圖。本文冷卻風(fēng)扇仿真模擬分析中湍流模型選擇 Realizable K-Epsilon Two-Layer （all Y+），該湍流模型要求壁面Y+值盡量在0～5之間。從圖9（a）可以看出，冷卻風(fēng)扇旋轉(zhuǎn)葉片Wall Y+值大部分區(qū)域滿足Wall Y+要求，表明本文冷卻風(fēng)扇的計(jì)算域網(wǎng)格滿足仿真分析要求。從圖9（b）和圖9（c）中可以看出，冷卻風(fēng)扇葉片表面速度分布較均勻，冷卻風(fēng)扇外圈因風(fēng)扇罩設(shè)計(jì)的原因?qū)е庐a(chǎn)生較嚴(yán)重的風(fēng)量回流，引起風(fēng)扇效率較低，與實(shí)驗(yàn)測(cè)試情況基本一致。后續(xù)研究中可以考慮對(duì)風(fēng)扇罩或風(fēng)扇進(jìn)行優(yōu)化設(shè)計(jì)，提升冷卻風(fēng)扇的運(yùn)作效率。

圖9 冷卻風(fēng)扇仿真計(jì)算結(jié)果

表1和圖10為不同測(cè)試風(fēng)壓下冷卻風(fēng)扇風(fēng)量仿真結(jié)果與實(shí)驗(yàn)測(cè)試結(jié)果的對(duì)比分析。從表1和圖10中可以看出，冷卻風(fēng)扇風(fēng)量的仿真結(jié)果與實(shí)驗(yàn)測(cè)試結(jié)果的變化趨勢(shì)基本一致，最大誤差約為5%，總體平均誤差約為3%，表明本文構(gòu)建的冷卻風(fēng)扇內(nèi)流場(chǎng)特性的計(jì)算方法是可行的，可以有效地預(yù)測(cè)冷卻風(fēng)扇在不同測(cè)試風(fēng)壓下的風(fēng)量值，為驗(yàn)證冷卻風(fēng)扇實(shí)驗(yàn)測(cè)試結(jié)果提供了分析支撐，同時(shí)為后續(xù)對(duì)冷卻風(fēng)扇葉片的優(yōu)化設(shè)計(jì)及風(fēng)扇罩的優(yōu)化設(shè)計(jì)提供了理論基礎(chǔ)。

表1 風(fēng)量計(jì)算結(jié)果

圖10 冷卻風(fēng)扇風(fēng)量仿真結(jié)果與實(shí)驗(yàn)測(cè)試結(jié)果對(duì)比圖

5 結(jié)束語(yǔ)

（1）通過網(wǎng)格無關(guān)性研究，得出該冷卻風(fēng)扇計(jì)算域模型的網(wǎng)格控制策略為：非旋轉(zhuǎn)域網(wǎng)格尺寸16 mm，過渡區(qū)域網(wǎng)格尺寸8 mm，旋轉(zhuǎn)域網(wǎng)格尺寸1 mm，葉尖加密尺寸0.25 mm，邊界層層數(shù)6層，邊界層增長(zhǎng)比1.5，風(fēng)扇邊界層厚度0.4 mm，其余邊界層厚度1 mm，入口球面半徑為風(fēng)扇外輪廓半徑的4.2倍。該網(wǎng)格控制策略下可獲得較好地網(wǎng)格質(zhì)量，且風(fēng)扇旋轉(zhuǎn)葉片Wall Y+值滿足湍流模型的Wall Y+范圍。

（2）Realizable K-Epsilon Two-Layer （all Y+）湍流模型下冷卻風(fēng)扇風(fēng)量仿真結(jié)果與實(shí)測(cè)風(fēng)量值間的誤差較小，氣體物性采用常密度設(shè)置和采用穩(wěn)態(tài)求解器，可降低計(jì)算成本，且仿真結(jié)果與實(shí)驗(yàn)更吻合。

（3）基于計(jì)算流體動(dòng)力學(xué)方法，采用STAR-CCM+軟件對(duì)冷卻風(fēng)扇內(nèi)流場(chǎng)進(jìn)行了流場(chǎng)特性研究，冷卻風(fēng)扇風(fēng)量仿真結(jié)果與風(fēng)扇實(shí)驗(yàn)測(cè)試結(jié)果的變化趨勢(shì)一致，兩者間的誤差在可接受的范圍內(nèi)，從而驗(yàn)證了CFD仿真分析的可靠性。基于計(jì)算流體動(dòng)力學(xué)理論的數(shù)值仿真分析方法可以較好地預(yù)測(cè)冷卻風(fēng)扇內(nèi)部流場(chǎng)特性及其風(fēng)量值，研究結(jié)果可為后續(xù)冷卻風(fēng)扇及風(fēng)扇罩的結(jié)構(gòu)設(shè)計(jì)及性能優(yōu)化提供理論基礎(chǔ)。