不同網(wǎng)格劃分方式的圓柱型舵桿繞流數(shù)值模擬結(jié)果比較

周上然

（海軍上海地區(qū)裝備修理監(jiān)修室，上海 200000）

船舶在航行時(shí)，流體經(jīng)過(guò)舵桿會(huì)產(chǎn)生渦。雷諾數(shù)在一定范圍時(shí)，流體在舵桿表面的后駐點(diǎn)附近脫落，脫落漩渦以周期性的交替方式離開(kāi)圓柱表面，形成渦泄。渦泄的存在使得舵桿表面周期性變化的升力和阻力增加，從而導(dǎo)致物體的振動(dòng)，影響了舵桿的穩(wěn)定性。本文將舵桿簡(jiǎn)化為二維圓柱模型，運(yùn)用兩種網(wǎng)格劃分方式模擬Re=200時(shí)二維靜止圓柱的繞流情況，對(duì)模擬出的升力與阻力系數(shù)進(jìn)行比較與分析。

1 數(shù)值建模

1.1 基本方程及邊界條件

對(duì)不可壓縮粘性流體，在直角坐標(biāo)系下，其運(yùn)動(dòng)規(guī)律可用N-S方程來(lái)描述，連續(xù)性方程和動(dòng)量方程分別為

式中：

ui、uj分別為 x、y 方向的速度；

ρ為空氣密度；

p為壓力；

ν為流體動(dòng)力粘度系數(shù)。

選取整個(gè)流場(chǎng)區(qū)域大小為30 D×55 D(D為立管直徑，D=0.02 m，ν=0.010 03 m/s)，其中，上游為15 D，下游為40 D，立管圓心距上下兩側(cè)各15 D，左側(cè)、上側(cè)與下側(cè)設(shè)為速度入口，流體自左向右流動(dòng)，右側(cè)為速度出口，圓柱壁面滿足無(wú)滑移固壁條件。本算例中雷諾數(shù)為200為層流流動(dòng)，因此，在粘性模型的選取上采用Laminar層流模型。在計(jì)算方法上采用SIMPLEC方法，由于流動(dòng)為非定常流動(dòng)，為模擬瞬態(tài)流動(dòng)，時(shí)間離散方式采用二階隱式，對(duì)壓力項(xiàng)離散格式采用二階精度。

1.2 網(wǎng)格劃分

流場(chǎng)計(jì)算區(qū)域選為55 D×30 D，圓柱距入口邊界為15 D，上下邊界距圓柱為15 D，D為物體垂直于來(lái)流方向平面上的特征尺寸，對(duì)圓柱一般取直徑，D=0.02 m，ν=0.010 03 m/s。

網(wǎng)格高度對(duì)于計(jì)算結(jié)果是否精確有著決定性的作用，對(duì)于層流模型，圓柱近壁處第一層網(wǎng)格高度需要滿足y+≈1，根據(jù)如下公式估算第一層網(wǎng)格控制點(diǎn)離開(kāi)壁面距離△y：

第一種方法對(duì)計(jì)算區(qū)域進(jìn)行了劃分，中心圓柱采用O型網(wǎng)格并進(jìn)行局部加密，右側(cè)采用漸變的結(jié)構(gòu)網(wǎng)格以提高計(jì)算效率，如圖1所示。第二種方法將全流場(chǎng)分為9個(gè)區(qū)域，采用四邊形結(jié)構(gòu)網(wǎng)格對(duì)中間圓柱周?chē)?D區(qū)域及尾流區(qū)進(jìn)行局部加密，如圖2所示。

圖1 網(wǎng)格劃分

圖2 網(wǎng)格劃分

1.3 計(jì)算結(jié)果及分析

圖3、4分別為第一種與第二種網(wǎng)格劃分（Re=200）時(shí)的渦街脫落圖。

圖3 渦街脫落圖

圖4 渦街脫落圖

兩種網(wǎng)格模擬的現(xiàn)象清晰，基本相同，可以認(rèn)定兩種方式均可以模擬光滑圓柱舵桿繞流。

2 模擬結(jié)果的比較與分析

2.1 升力系數(shù)與阻力系數(shù)的分析

圖3、圖4的現(xiàn)象清晰的模擬了光滑立管尾流場(chǎng)，因此，可認(rèn)為兩種劃分方式均可行。圖5、圖6分別給出了兩種不同劃分格式的圓柱阻力系數(shù)與升力系數(shù)的時(shí)程曲線。從圖上可以看出圓柱尾流形成了卡門(mén)渦街，流動(dòng)在150 s后趨于穩(wěn)定，當(dāng)位于圓柱上側(cè)的旋渦脫落后，圓柱上方壓力增加，此時(shí)圓柱下方旋渦還未脫落，壓力較小，上下兩側(cè)產(chǎn)生負(fù)的壓力差，反之當(dāng)位于圓柱下側(cè)的旋渦脫落時(shí)，產(chǎn)生正的壓力差，因此，升力系數(shù)以y=0為平衡位置進(jìn)行周期性震蕩，升力系數(shù)的均值為0。升力系數(shù)與阻力系數(shù)曲線在渦穩(wěn)定脫落后，均呈現(xiàn)周期性變化[1]。兩種劃分方式計(jì)算出的結(jié)果：

（1）升力系數(shù)幅值為0.51，阻力系數(shù)在渦脫達(dá)到穩(wěn)定后的均值為1.23。

（2）升力系數(shù)幅值為0.53，阻力系數(shù)在渦脫達(dá)到穩(wěn)定后的均值為1.32。

從升力系數(shù)和阻力系數(shù)的周期性變化中，可以清楚的反映出圓柱上旋渦脫落的特性。均勻來(lái)流經(jīng)過(guò)圓柱表面時(shí)，流動(dòng)發(fā)生變化，產(chǎn)生邊界層分離后，在逆壓梯度的作用下，于圓柱上下兩個(gè)表面附近產(chǎn)生旋渦，旋渦在來(lái)流的作用下隨著流體向后運(yùn)動(dòng)，然后脫落。由于圓柱上下兩側(cè)的旋渦變化對(duì)流動(dòng)方向上的壓差影響是一樣的，所以來(lái)流壓力差的變化，較垂直來(lái)流方向的壓力差變化快，變化速度是它的兩倍，這就很好的解釋了阻力系數(shù)的振蕩周期為升力系數(shù)振蕩周期的二分之一。

圖5 第一種劃分方式

圖6 第二種劃分方式

對(duì)升力系數(shù)時(shí)程曲線做傅里葉變換，從時(shí)域轉(zhuǎn)換到頻域，第一種劃分方式計(jì)算出的主頻率為0.095 Hz，如圖7所示，第二種劃分方式計(jì)算出的主頻率為0.102 Hz，如圖8所示。

圖7 第一種劃分方式升力系數(shù)時(shí)程譜

圖8 第二種劃分方式升力系數(shù)時(shí)程譜

2.2 斯特勞哈數(shù)(St)

St數(shù)是描述圓柱繞流的一個(gè)重要參數(shù)，其將邊界層分離及流動(dòng)的不穩(wěn)定性這些微觀的隨機(jī)特性，與表現(xiàn)相對(duì)穩(wěn)定性的渦泄頻率聯(lián)系在一起[2]，其定義為

式中：

St為Strouhal數(shù)，表示旋渦的脫離情況；

f為圓柱渦泄頻率；

U為未受干擾的自由來(lái)流速度。

一般而言，在絕大多數(shù)的Re數(shù)區(qū)域內(nèi)，St數(shù)相對(duì)穩(wěn)定，這意味著渦泄有序。第一種劃分方式計(jì)算出Strouhal數(shù)為0.19，第二種劃分方式計(jì)算出Strouhal數(shù)為0.2。

2.3 結(jié)果比較

表1中給出了國(guó)內(nèi)外的研究學(xué)者對(duì)雷諾數(shù)為200時(shí)的圓柱繞流數(shù)值計(jì)算[3]與實(shí)驗(yàn)結(jié)果，通過(guò)比較，第二種劃分方式的計(jì)算出的Cl與Cd的值與實(shí)驗(yàn)結(jié)果相接近，可以認(rèn)為第二種劃分方式計(jì)算精度較高。

表1 Re=200時(shí)圓柱繞流數(shù)值計(jì)算結(jié)果與文獻(xiàn)數(shù)據(jù)比較

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

通過(guò)以上分析，我們得出以下結(jié)論：

（1）對(duì)Re=200的低雷諾數(shù)圓柱舵桿繞流，采用層流模型并結(jié)合合理的網(wǎng)格劃分，可以很好地模擬舵桿渦流運(yùn)動(dòng)，并通過(guò)對(duì)升力系數(shù)、阻力系數(shù)的分析，得到渦的脫落頻率。

（2）雖然兩種劃分方式都可以模擬出圓柱舵桿繞流運(yùn)動(dòng)，但是第二種方法計(jì)算出的升力系數(shù)、阻力系數(shù)更加接近實(shí)驗(yàn)結(jié)果，因?yàn)椋湓趤?lái)流方向與垂直來(lái)流方向均加密了網(wǎng)格，較好的覆蓋了圓柱舵桿的分離區(qū)域。

（3）本文只對(duì)靜止的二維圓柱模型作了數(shù)值計(jì)算，未考慮圓柱的運(yùn)動(dòng)對(duì)流場(chǎng)產(chǎn)生的影響與三維模型的尺度效應(yīng)。

[1]李 田，張繼業(yè)，等.二維彈性圓柱渦致振動(dòng)尾渦模態(tài)[J].空氣動(dòng)力學(xué)學(xué)報(bào)，2010，(06):25-26.

[2]顧 罡.二維圓柱繞流、雙圓柱繞流問(wèn)題和三維垂蕩板運(yùn)動(dòng)的數(shù)值模擬[J].上海交通大學(xué)，2007，23-24.

[3]苑健康，黃維平.二維立管Helical Strakes繞流場(chǎng)的ANSYSCFD 分析[J].船海工程，2010，(04):26-27.