李邦華，黎 峰，晉文菊，劉向冬，蔣曙暉

(上海船舶研究設(shè)計(jì)院，上海 201203)

0 引 言

當(dāng)流體在近壁面流動(dòng)時(shí)，近壁區(qū)域各物理變量的變化梯度非常大，因此該區(qū)域在采用CFD 進(jìn)行數(shù)值計(jì)算時(shí)需要著重考慮。對(duì)于近壁區(qū)域的處理目前主要有2 種方法：一種是壁面函數(shù)法，該方法無需對(duì)過渡層和粘性底層進(jìn)行直接求解，而是通過半經(jīng)驗(yàn)公式將湍流核心區(qū)域的求解變量與壁面上的物理量聯(lián)系起來；另一種方法是近壁模型法，使用該方法時(shí)需要在壁面附近劃分足夠細(xì)的網(wǎng)格，使之可以求解粘性影響的區(qū)域[1]。

對(duì)于k-ε 兩方程高雷諾數(shù)湍流模型與壁面函數(shù)相結(jié)合的模式，邊界層網(wǎng)格尺度對(duì)計(jì)算結(jié)果有較大影響，即需要保證y+值在一個(gè)合理的范圍[2]。資丹等[3]提出3 種具有不同邊界層網(wǎng)格尺度的網(wǎng)格劃分方案，分析了不同初始y+值對(duì)泵站流場數(shù)值計(jì)算結(jié)果的影響。杜云龍等[4]確定出船舶阻力數(shù)值計(jì)算中對(duì)結(jié)果影響最大的因素為邊界層網(wǎng)格尺度。Jialun Liu 等[5]系統(tǒng)研究了網(wǎng)格類型、計(jì)算域形狀、網(wǎng)格大小以及網(wǎng)格增長率對(duì)船舶舵系水動(dòng)力計(jì)算結(jié)果的影響。Nguyen Van Trieu[6]采用SST k-w 壁面模型法研究魚尾舵隨邊寬度對(duì)水動(dòng)力性能的影響。賴晨光等[7]研究了不同湍流模型在匹配不同邊界層網(wǎng)格策略時(shí)對(duì)汽車外流場數(shù)值計(jì)算精確性的影響。

本文所進(jìn)行的高效舵水動(dòng)力數(shù)值計(jì)算屬于高雷諾數(shù)粘性繞流問題，采用壁面函數(shù)法就可以有效地模擬流體的流動(dòng)。為了研究邊界層網(wǎng)格尺度在高效舵水動(dòng)力計(jì)算中的影響，在網(wǎng)格數(shù)量無關(guān)性檢查的基礎(chǔ)上，針對(duì)幾種典型的初始壁面y+值，應(yīng)用Standard k-ε 湍流模型對(duì)一種新型隨邊扭曲高效舵進(jìn)行水動(dòng)力計(jì)算，并將模擬結(jié)果與模型試驗(yàn)結(jié)果進(jìn)行對(duì)比分析，從而為高效舵水動(dòng)力計(jì)算提供參考。

1 模型試驗(yàn)

舵的敞水模型試驗(yàn)在上海船舶運(yùn)輸科學(xué)研究所的拖曳水池中進(jìn)行。水池長192 m，深4.2 m，寬10 m，設(shè)置在水池上的拖車系統(tǒng)運(yùn)行速度為0.01～10 m/s。采用六分力傳感器測量作用在舵上的升力FL、阻力Fd和扭矩M。

試驗(yàn)工況和縮尺比如表1 所示。本試驗(yàn)需測量在不同舵角下舵的受力情況，力的測量從0°～30°范圍每隔5°進(jìn)行一次。試驗(yàn)時(shí)，首先需調(diào)整舵角至相應(yīng)角度，再將拖車開至工況要求的速度，最后拖車保持勻速行駛至試驗(yàn)結(jié)束。舵的尾流場如圖1 所示。

2 控制方程與湍流模型

2.1 控制方程

本文在數(shù)值計(jì)算時(shí)假定流體為密度不變，不可壓縮的粘性流，其連續(xù)性方程和動(dòng)量方程可以簡化如下：

圖 1 舵的尾流場Fig. 1 The wake field of the rudder

式中：u，v，w 分別為在x，y，z 方向上的速度分量；g 為重力加速度；ρ 為流體密度；為動(dòng)力粘度；流體內(nèi)部壓強(qiáng)為p；Sx，Sy，Sz分別為x，y，z 方向上的動(dòng)量源項(xiàng)。

2.2 湍流模型

Standard k-ε 湍流模型[8]是較為常用的兩方程湍流模型，其輸運(yùn)方程如下：

其中：Gk為平均速度引起的湍動(dòng)能產(chǎn)生項(xiàng)；Gb為浮力引起的湍動(dòng)能產(chǎn)生項(xiàng)；渦黏系數(shù)為ut，Ym表示湍動(dòng)耗散率；C1ε，C2ε，C3ε，σε和σk為經(jīng)驗(yàn)常數(shù)（取值參考文獻(xiàn)[8]）。

3 計(jì)算域的建立

舵的幾何模型與模型試驗(yàn)保持一致，模型的建立在Rhinoceros 中完成。舵的外流場為一個(gè)單獨(dú)的長方體計(jì)算域，如圖2 所示。舵前端距速度入口邊界4 倍弦長，來流速度與模型試驗(yàn)拖車速度一致；舵隨邊距壓力出口邊界7 倍弦長；計(jì)算域?qū)挾确较蛉? 倍弦長，左右側(cè)面設(shè)定為壁面邊界；舵頂和舵底分別距流場底部和頂部4 倍舵高，為壁面邊界；舵表面設(shè)為無滑移壁面邊界。模型試驗(yàn)時(shí)，測得的水溫為19.2 C°，因此數(shù)值模擬時(shí)水的動(dòng)力粘度為0.001 029 9 Pa·s，密度取為998.26 kg/m3。

圖 2 計(jì)算域及邊界面Fig. 2 Computational domain and the boundary condition

4 數(shù)值計(jì)算分析

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

網(wǎng)格的劃分情況對(duì)流體的運(yùn)動(dòng)以及數(shù)值計(jì)算的精度和收斂情況具有重要的影響[9–10]，因此在探討邊界層網(wǎng)格尺度對(duì)高效舵水動(dòng)力數(shù)值計(jì)算影響前，對(duì)邊界層外的網(wǎng)格進(jìn)行無關(guān)性驗(yàn)證。本次驗(yàn)證對(duì)計(jì)算域采用4 種不同的網(wǎng)格尺度，分別對(duì)應(yīng)4 種不同數(shù)量的網(wǎng)格mesh1，mesh2，mesh3，mesh4 和mesh5。計(jì)算時(shí)采用Standard k-ε 湍流模型和相同的邊界層網(wǎng)格尺度。網(wǎng)格無關(guān)性驗(yàn)證結(jié)果如表2 和圖3 所示。

表 2 四種不同網(wǎng)格數(shù)量下的升力系數(shù)對(duì)比Tab. 2 Comparison of lift coefficient with different mesh scale

圖 3 不同網(wǎng)格數(shù)量下的升力系數(shù)對(duì)比Fig. 3 Comparison of lift coefficient with different mesh scale

從表2 和圖3 可以看出，隨著網(wǎng)格數(shù)量的增加數(shù)值計(jì)算結(jié)果變化不大，同時(shí)由mesh2 和mesh4 對(duì)比可以看出，當(dāng)網(wǎng)格數(shù)量達(dá)到152 萬后，增加網(wǎng)格數(shù)量對(duì)計(jì)算結(jié)果沒有明顯的影響，不過較為細(xì)密的網(wǎng)格對(duì)舵周圍流場細(xì)節(jié)的捕捉更為精確，因此最終選擇mesh3為后續(xù)數(shù)值計(jì)算的網(wǎng)格劃分方案。

4.2 邊界層網(wǎng)格尺度影響

舵在流場中運(yùn)動(dòng)時(shí)，由于水具有粘性，流體在近壁區(qū)域的流動(dòng)大致可以分為內(nèi)、中、外3 層：貼近壁面的一層稱為粘性底層，該層的流動(dòng)幾乎是層流；過渡層處于粘性底層的外面，流體的流動(dòng)狀態(tài)介于層流與湍流之間；最外層是完全湍流層也稱為對(duì)數(shù)層，湍流處于充分發(fā)展?fàn)顟B(tài)。為了合理描述邊界層內(nèi)的流動(dòng)，引入體現(xiàn)邊界層網(wǎng)格尺度的無量綱參數(shù)y+，其表達(dá)式如公式4 所示[11]。通常認(rèn)為，當(dāng)y+＜5 的區(qū)域?yàn)檎承缘讓樱?＜y+＜30 的區(qū)域?yàn)檫^渡層，y+＞30 的區(qū)域?yàn)閷?duì)數(shù)層。

其中：L 為特征長度，Re 為雷諾數(shù)，△yp為第1 層網(wǎng)格高度。

為了研究邊界層網(wǎng)格尺度對(duì)數(shù)值模擬結(jié)果產(chǎn)生的影響，選用典型的初始壁面y+值3，15，30，60 和150 對(duì)15°攻角下的高效舵進(jìn)行水動(dòng)力計(jì)算，其中邊界層內(nèi)的網(wǎng)格增長率取1.2。不同初始壁面y+值下數(shù)值計(jì)算結(jié)果與模型試驗(yàn)結(jié)果對(duì)比如表3 所示。

表 3 不同初始壁面y+值下數(shù)值計(jì)算結(jié)果與模型試驗(yàn)結(jié)果對(duì)比Tab. 3 Comparison of model test with numerical calculations at different y+

從表3 和圖4 的對(duì)比可以看出，隨著初始壁面y+值的增加，升力系數(shù)與模型試驗(yàn)誤差從–0.807%變化到2.033%，阻力系數(shù)從4.622%變化到–4.507%，從而說明初始壁面y+值對(duì)高效舵水動(dòng)力計(jì)算精度有著較為明顯的影響；阻力系數(shù)對(duì)初始壁面y+值的敏感性較升力系數(shù)更為明顯，這是由于初始壁面y+值直接反應(yīng)了邊界層內(nèi)網(wǎng)格尺度，邊界層內(nèi)的網(wǎng)格尺度直接影響舵近壁區(qū)速度梯度的變化，從而影響作用在舵上的阻力，但由于壁面函數(shù)的應(yīng)用第1 層邊界層厚度不宜過薄，即y+不宜過小。從圖4 還可以看出，隨著y+值的增加，舵尾靜流區(qū)的范圍逐漸減小，舵周圍流場速度分布也略有不同，這也能說明初始壁面y+值對(duì)舵水動(dòng)力計(jì)算會(huì)產(chǎn)生影響。綜合來看當(dāng)初始壁面y+值取30 時(shí)，數(shù)值計(jì)算結(jié)果與模型試驗(yàn)結(jié)果誤差最小。

4.3 不同攻角下的數(shù)值計(jì)算結(jié)果

圖 5 數(shù)值模擬與模型試驗(yàn)結(jié)果對(duì)比Fig. 5 Comparison of numerical calculations with model test

基于前文的網(wǎng)格劃分方案和壁面初始y+值，對(duì)隨邊扭曲高效舵進(jìn)行攻角從0°～30°（間隔為5°）下數(shù)值模擬，并將模型試驗(yàn)與模擬結(jié)果進(jìn)行了對(duì)比，如圖5所示?？梢钥闯?，數(shù)值模擬與模型試驗(yàn)結(jié)果吻合度較高，尤其在失速角以前。

5 結(jié) 語

為了研究邊界層網(wǎng)格尺度在高效舵水動(dòng)力計(jì)算中的影響，在網(wǎng)格數(shù)量無關(guān)性檢查的基礎(chǔ)上，針對(duì)幾種典型的y+值，應(yīng)用Standard k-ε 湍流模型對(duì)一種新型隨邊扭曲高效舵進(jìn)行水動(dòng)力計(jì)算，并與模型試驗(yàn)結(jié)果進(jìn)行對(duì)比。

可以看出，初始壁面y+值對(duì)高效舵水動(dòng)力計(jì)算精度有明顯的影響；阻力系數(shù)對(duì)初始壁面y+值的敏感性較升力系數(shù)更為顯著；當(dāng)初始壁面y+值取30 時(shí)，數(shù)值計(jì)算結(jié)果與模型試驗(yàn)結(jié)果誤差最小，因此建議進(jìn)行高效舵水動(dòng)力數(shù)值計(jì)算時(shí)初始壁面y+值取在30 左右。