零航速減搖鰭升力特性數(shù)值研究

周　衛(wèi)，張吉萍，許頌捷

（1.浙江海洋學(xué)院船舶與海洋工程學(xué)院，浙江舟山　316022；2.浙江省海洋漁業(yè)船舶交易服務(wù)中心，浙江杭州　310007；3.浙江省近海海洋工程技術(shù)重點(diǎn)實(shí)驗(yàn)室，浙江舟山　316022；4.上海羽翼船舶設(shè)備有限公司，上海　201600）

零航速減搖鰭升力特性數(shù)值研究

周衛(wèi)1，2，張吉萍1，3，許頌捷4

（1.浙江海洋學(xué)院船舶與海洋工程學(xué)院，浙江舟山316022；2.浙江省海洋漁業(yè)船舶交易服務(wù)中心，浙江杭州310007；3.浙江省近海海洋工程技術(shù)重點(diǎn)實(shí)驗(yàn)室，浙江舟山316022；4.上海羽翼船舶設(shè)備有限公司，上海201600）

為了探究零航速減搖鰭周期擺動(dòng)時(shí)產(chǎn)生的升力與轉(zhuǎn)動(dòng)頻率、轉(zhuǎn)角幅度、轉(zhuǎn)軸位置以及來(lái)流速度幾方面的關(guān)系，本文運(yùn)用滑移網(wǎng)格技術(shù)對(duì)二維NACA0012翼型的周期擺動(dòng)進(jìn)行數(shù)值預(yù)報(bào)。研究表明：Transition SST湍流模型較其他常見(jiàn)湍流模型更適合分析大攻角周期擺動(dòng)問(wèn)題；對(duì)于轉(zhuǎn)動(dòng)方式滿(mǎn)足余弦函數(shù)的減搖鰭，產(chǎn)生的升力隨轉(zhuǎn)角幅度、轉(zhuǎn)動(dòng)頻率增加而增大，來(lái)流速度對(duì)升力的影響大于轉(zhuǎn)軸位置。論文研究對(duì)零航速減搖鰭選型設(shè)計(jì)具有一定指導(dǎo)意義。

零航速；減搖鰭；升力特性；Transition SST；NACA0012

隨著人們對(duì)低航速下的舒適性要求的提高，零航速減搖概念自提出便受到國(guó)內(nèi)外機(jī)構(gòu)及學(xué)者的關(guān)注，其工作原理是快速擺動(dòng)鰭體，產(chǎn)生抵抗船舶橫搖力矩，實(shí)現(xiàn)低航速或停泊狀態(tài)的減搖。目前關(guān)于零航速減搖鰭的研究主要以試驗(yàn)研究為主，美國(guó)及荷蘭相關(guān)機(jī)構(gòu)在1998年開(kāi)展了探索性試驗(yàn)研究［1-2］，昆騰公司利用試驗(yàn)方法開(kāi)展全航速減搖鰭的初步設(shè)計(jì)［3］。在國(guó)內(nèi)，金鴻章等［4］根據(jù)實(shí)驗(yàn)數(shù)據(jù)及勢(shì)流理論構(gòu)建零航速減搖鰭的力學(xué)模型，表明升力受旋轉(zhuǎn)角速度、角加速度等因素的影響。

零航速減搖鰭的周期性擺動(dòng)實(shí)際是以旋渦產(chǎn)生、耗散為主導(dǎo)的湍流運(yùn)動(dòng)，對(duì)于這種復(fù)雜的湍流運(yùn)動(dòng)，試驗(yàn)研究受試驗(yàn)器材和設(shè)備的限制，往往僅能提供有限的積分變量，而CFD數(shù)值方法能提供更完善的流場(chǎng)信息。同時(shí)伴隨湍流模型的完善，及計(jì)算精度與穩(wěn)定性的提高，運(yùn)用數(shù)值方法進(jìn)行零航速減搖鰭的研究具有較明顯的優(yōu)勢(shì)。

本文以二維NACA0012翼型作為研究對(duì)象，運(yùn)用Fluent軟件自帶的UDF函數(shù)編譯各種運(yùn)動(dòng)方式，并結(jié)合滑移網(wǎng)格技術(shù)進(jìn)行數(shù)值預(yù)報(bào)。針對(duì)試驗(yàn)工況，對(duì)各類(lèi)湍流模型、網(wǎng)格大小等相關(guān)設(shè)置的影響進(jìn)行分析，表明數(shù)值方法能較好的預(yù)報(bào)大攻角動(dòng)態(tài)失速問(wèn)題。其次針對(duì)減搖鰭運(yùn)動(dòng)特征，分析了轉(zhuǎn)角幅度、轉(zhuǎn)動(dòng)周期、轉(zhuǎn)軸位置、來(lái)流速度幾方面對(duì)于升力的影響。

1　數(shù)值方法

1.1湍流模型

直接計(jì)算法和雷諾平均法（RANS）是進(jìn)行數(shù)值分析的兩種主要方法，RANS由于存在表示湍流脈動(dòng)對(duì)平均流動(dòng)影響的雷諾應(yīng)力項(xiàng)τij，因此需要添加湍流模型進(jìn)行封閉。常用的k-ω和k-ε湍流模型，由于存在較大耗散會(huì)制約旋渦的產(chǎn)生，無(wú)法準(zhǔn)確的預(yù)報(bào)旋渦的產(chǎn)生的位置、旋渦擴(kuò)散的過(guò)程［5］。因此，本文主要采用Transition SST模型進(jìn)行分析減搖鰭的周期擺動(dòng)問(wèn)題。

Transition SST湍流模型是基于k-ω SST模型，同時(shí)考慮湍流轉(zhuǎn)捩問(wèn)題，因此運(yùn)輸方程還需構(gòu)造關(guān)于動(dòng)能厚度雷諾數(shù)Reθ和間歇因子γ的運(yùn)輸方程，并通過(guò)控制方程與SST k-ω模型結(jié)合起來(lái)。

式（1）、（2）為SST k-ω湍流模型，式（3）、（4）為間歇因子γ及轉(zhuǎn)捩雷諾數(shù)的運(yùn)輸方程；式（1）、（2）中：ρ為流體密度，μj對(duì)應(yīng)速度矢量，xj對(duì)應(yīng)位置矢量，Gk和Gω為關(guān)于k和ω的產(chǎn)生項(xiàng)，Гω、Гk為有效擴(kuò)散項(xiàng)；Yω、Yk為各自的湍流耗散項(xiàng)，而Sω、Sk則表示為可自定義的源項(xiàng)。式（3）中Pγ1、Eγ1為過(guò)渡源項(xiàng)，Pγ2、Eγ2為破壞源項(xiàng)。

轉(zhuǎn)捩模型與SST湍流模型結(jié)合的控制方程為：

G%k與Y*k為SST模型中的初始的生成項(xiàng)與耗散項(xiàng)。

1.2研究對(duì)象

目前，關(guān)于零航速減搖鰭的試驗(yàn)數(shù)據(jù)有限，因此文本選用LEE等［6］針對(duì)大攻角失速研究采用的NACA0012翼型作為研究對(duì)象。選用該試驗(yàn)數(shù)據(jù)進(jìn)行驗(yàn)證分析主要有兩方面原因：首先是該實(shí)驗(yàn)涉及大攻角翼型的雷諾數(shù)范圍與減搖鰭擺動(dòng)過(guò)程的雷諾數(shù)范圍相近，同時(shí)空氣與水的流動(dòng)特性相近。其次，該試驗(yàn)得到的各項(xiàng)系數(shù)已剔除三維流動(dòng)對(duì)是實(shí)驗(yàn)結(jié)果的影響，因此試驗(yàn)數(shù)據(jù)適用于二維的分析。NACA0012翼型的擺動(dòng)原理具體設(shè)置可參見(jiàn)文獻(xiàn)［6］。表1為數(shù)值計(jì)算涉及工況的主要參數(shù)。

表1　工況參數(shù)Tab.1　Parameters of Different Loadcases

1.3網(wǎng)格及計(jì)算設(shè)置

文中所計(jì)算的兩種工況采用相同的網(wǎng)格進(jìn)行計(jì)算，網(wǎng)格利用ICEM軟件的O型網(wǎng)格形式進(jìn)行劃分，翼型表面初始網(wǎng)格高度為0.02mm，滿(mǎn)足y+值小于1要求，網(wǎng)格沿法線(xiàn)方向的高度以1.2比率增長(zhǎng)，確保能捕獲近壁面附近的小尺度脈動(dòng)同時(shí)降低網(wǎng)格數(shù)量，最終確定的初始網(wǎng)格數(shù)量為12萬(wàn)。在此基礎(chǔ)上，弦長(zhǎng)方向上的網(wǎng)格加密1.5倍，分別生成18萬(wàn)和27萬(wàn)網(wǎng)格。

圖1　NACA翼型表面網(wǎng)格Fig.1　Surface Mesh of NACA Airfoils

圖2　計(jì)算域及邊界條件Fig.2　Calculation Domain and boundary conditions

數(shù)值計(jì)算采用的邊界條件與試驗(yàn)條件相同，具體設(shè)置如下：入口位置：翼前18倍弦長(zhǎng)，設(shè)置為速度入口；遠(yuǎn)場(chǎng)條件：距翼型18倍弦長(zhǎng)，設(shè)置為速度入口；出口位置：翼后24倍弦長(zhǎng)位置，設(shè)置為壓力出口；Interface邊界：距翼型10倍弦長(zhǎng)。

本文涉及的計(jì)算均采用雙精度基于壓力法求解器進(jìn)行求解，對(duì)流項(xiàng)采用二階迎風(fēng)差分格式，壓力速度差分方程組的偶合求解采用SIMPLE格式，采用中心差分格式進(jìn)行離散。對(duì)于DDES湍流模型，采用有界中心差分格式。運(yùn)用滑移網(wǎng)格技術(shù)以及UDF函數(shù)實(shí)現(xiàn)不同運(yùn)動(dòng)方式的控制。

2　計(jì)算結(jié)果與分析

2.1驗(yàn)證工況分析

圖3為時(shí)間步長(zhǎng)是0.000 2 s下，不同湍流模型計(jì)算得到升力系數(shù)。對(duì)于所有湍流模型對(duì)于從-5°到25°的上揚(yáng)階段準(zhǔn)確性均優(yōu)于25°到-5°的下擺階段，無(wú)論是升力變化趨勢(shì)還是與試驗(yàn)結(jié)果的差值。

上擺階段：角度滿(mǎn)足-5°≤α≤25°范圍時(shí)：RNG k-ε、DDES、Transition SST三種湍流模型均做出較好的預(yù)報(bào)，在20°≤α≤25°范圍，三種湍流模型表現(xiàn)較大差異，DDES模型預(yù)報(bào)的失速角度為21°，且升力系數(shù)峰值較實(shí)驗(yàn)值偏小15.6%；RNG k-ε模型則升力系數(shù)隨轉(zhuǎn)角增大而保持震蕩增加，并未出現(xiàn)明顯的失速現(xiàn)象，這是由于RNG k-ε模型無(wú)法模擬前緣渦分離現(xiàn)象。Transition SST模型預(yù)測(cè)的失速角度為24.5°，25°時(shí)的升力系數(shù)下降約0.9，表明Transition SST模型過(guò)度預(yù)報(bào)分離現(xiàn)象，將直接影響下擺過(guò)程的準(zhǔn)確性。

下擺階段：采用加強(qiáng)避免函數(shù)處理的RNG k-ε得到的升力系數(shù)變化較為平穩(wěn)，與試驗(yàn)結(jié)果吻合較好；Transition SST模型和DDES模型在15°≤α≤25°范圍，升力系數(shù)波蕩較大，這是由于過(guò)度預(yù)報(bào)分離現(xiàn)象，造成的二次分離預(yù)報(bào)不準(zhǔn)確。在15°≤α≤25°，Transition SST模型預(yù)報(bào)的升力系數(shù)的震蕩現(xiàn)象減小明顯，趨勢(shì)與實(shí)驗(yàn)結(jié)果吻合。

圖4為12萬(wàn)、18萬(wàn)、27萬(wàn)網(wǎng)格采用Transition SST模型、時(shí)間步長(zhǎng)采用0.0002 s時(shí)計(jì)算得到的升力系數(shù)。上揚(yáng)階段：三種網(wǎng)格得到的升力系數(shù)幾乎貼合，在20°至25°過(guò)程，12萬(wàn)網(wǎng)格與其他網(wǎng)格得到升力系數(shù)數(shù)值差異較大；下擺階段：三種網(wǎng)格得到升力系數(shù)變化趨勢(shì)整體保持一致，27萬(wàn)網(wǎng)格與試驗(yàn)結(jié)果的數(shù)值差異最小。整體而言，18萬(wàn)網(wǎng)格策略為較優(yōu)的網(wǎng)格策略，在后續(xù)的計(jì)算分析中采用18萬(wàn)網(wǎng)格進(jìn)行研究分析。

圖3　不同湍流模型的升力系數(shù)Fig.3　Lift Parameters of Different Turbulence Models

圖4　不同網(wǎng)格劃分的升力系數(shù)Fig.4　Lift Parameters of Different Mesh

2.2零航速數(shù)值分析

零航速減搖鰭的運(yùn)動(dòng)方式主要有兩種，一是采用加速-勻速-減速的運(yùn)動(dòng)方式［7］，二是滿(mǎn)足余弦運(yùn)動(dòng)規(guī)律，文中采用第二種運(yùn)動(dòng)方式進(jìn)行分析，并對(duì)轉(zhuǎn)角幅度、轉(zhuǎn)動(dòng)周期對(duì)升力影響進(jìn)行分析。

圖5　不同轉(zhuǎn)角幅度的升力變化曲線(xiàn)Fig.5　Lift Curve of different Angle Ranges

圖6　不同轉(zhuǎn)角周期的升力變化曲線(xiàn)Fig.6　Lift Curve of different Angle Periods

圖5為零航速減搖鰭的轉(zhuǎn)動(dòng)周期為1 s，轉(zhuǎn)角幅度分別為35°、40°、45°、50°時(shí)的升力變化曲線(xiàn)，由圖5可以看到，四種轉(zhuǎn)角幅度的升力曲線(xiàn)在上揚(yáng)或下擺階段升力峰值均經(jīng)歷先增大后減小的過(guò)程。下擺階段除轉(zhuǎn)角幅度為50°時(shí)升力峰值位于11°附近，其他轉(zhuǎn)角幅度的升力峰值出現(xiàn)在0°（最大角速度位置）附近；上揚(yáng)階段所有轉(zhuǎn)角幅度的升力峰值均出現(xiàn)在12°附近，即超前于最大角速度所在位置，這是由于附連水質(zhì)量力的影響。此外，轉(zhuǎn)角幅度增加5°，下擺階段升力峰值平均增幅約35%，上揚(yáng)階段升力峰值平均增幅約22%。

圖6對(duì)應(yīng)的是轉(zhuǎn)角幅度為40°，轉(zhuǎn)動(dòng)周期分別為0.5 s、1 s、2 s時(shí)的升力變化曲線(xiàn)。通過(guò)比較可以發(fā)現(xiàn)，三種周期下的升力峰值出現(xiàn)的位置基本一致，升力峰值與轉(zhuǎn)動(dòng)周期的平方成反比，即升力峰值與轉(zhuǎn)動(dòng)角速度平方成正比。

圖7　不同轉(zhuǎn)軸位置的升力變化曲線(xiàn)Fig.7　Lift Curve of different Rotor Positions

圖8　不同流速下的升力變化曲線(xiàn)Fig.8　Lift Curve of different Velocitys

圖7為四種轉(zhuǎn)軸位置的升力變化曲線(xiàn)，鰭軸位置向前移動(dòng)，升力變化趨勢(shì)基本保持不變，上擺過(guò)程的升力峰值增大，下擺過(guò)程的升力峰值減小。由圖8可以看到，流速達(dá)到0.2m/s時(shí)的升力變化趨勢(shì)和升力峰值均與航速低于0.1m/s時(shí)發(fā)生明顯變化，上擺和下擺過(guò)程升力變化趨勢(shì)相似、升力峰值相近，表明零、低航速時(shí)減搖鰭的升力變化受航速影響較大。

3　結(jié)論

本文利用數(shù)值方法對(duì)零航速減搖鰭上下擺動(dòng)過(guò)程中涉及的大攻角動(dòng)態(tài)失速問(wèn)題進(jìn)行研究，通過(guò)上述的分析可以得到以下結(jié)論：

（1）對(duì)于大攻角動(dòng)態(tài)失速問(wèn)題，Transition SST湍流模型的預(yù)報(bào)結(jié)果優(yōu)于常用的RNG k-ε湍流模型，DDES湍流模型存在過(guò)度預(yù)報(bào)旋渦分離現(xiàn)象。

（2）運(yùn)動(dòng)規(guī)律滿(mǎn)足α=α1cos（2πf·t）的零航速減搖鰭，在轉(zhuǎn)動(dòng)過(guò)程中減搖鰭產(chǎn)生的升力與轉(zhuǎn)角幅度成正比、同時(shí)與轉(zhuǎn)動(dòng)周期的平方成反比。改變鰭軸位置對(duì)于升力的變化影響不大，來(lái)流速度對(duì)于減搖鰭的性能影響明顯。

（3）對(duì)于零航速減搖鰭設(shè)計(jì)初期，可利用數(shù)值方法輔助開(kāi)展翼型評(píng)估、轉(zhuǎn)軸位置分析、轉(zhuǎn)動(dòng)形式優(yōu)化等方面工作。

［1］DALLINGA R P.Roll stabilization of motor yacht：Use of fin stabilizers in anchored conditions［Z］.Project 99，Amsterdam，1999.

［2］DALLINGA R P.Roll stabilization of motor yacht：Hydrodynamic aspects of the comparison of anti-roll tanks and fins［Z］. Project 2002，Amsterdam，2002.

［3］OOMS J.The use of roll stabilizer fins at zero speed［Z］.Project 2002，Amsterdam，2002.

［4］金鴻章，張曉飛.零航速減搖鰭升力模型研究［J］.海洋工程，2007，25（3）：83-87.

［5］ANSYS.ANSYS FLUENT 12.0 Theory Guide［Z］.Ansys Inc，USA，2009.

［6］LEE T，GERONTAKOS P.Investigation of flow over an oscillating airfoil［J］.Journal of Fluid Mechanics，2004，512：313-341.

［7］宋吉廣.基于升力反饋的全航速減搖鰭研究［D］.哈爾濱：哈爾濱工程大學(xué)，2010.

Research on Lift Characteristics of Zero-speed Fin Stabilizer

ZHOU Wei1，2，ZHANG Ji-ping1，3，XU Song-jie4
（1.School of Naval Architecture and Ocean Engineering，Zhejiang Ocean University，Zhoushan316022；2.Zhejiang Province Marine Fishing Vessel Trading Service Center，Hangzhou310007；3.Key Laboratory of Offshore Engineering Technology of Zhejiang Provine，Zhoushan316022；4.Yoyeah Marine，Shanghai201600，China）

To explore the influence from rotational frequency，angle range，rotor position and flow velocity to the cyclic swing lift of Zero-speed Fin Stabilizer，the numerical forecast to the cyclic swing of two-dimensional NACA0012 airfoils is carried out with the sliding mesh technique.The following conclusions can be obtained：Transition SST turbulence model is more appropriate to analysis the cyclic swing at high angles of attack that other turbulence models；to the fin stabilzer with the rotating mode meeting the cosine function，the lift increased with the increasing of angle range and rotational frequency，and the influence from the flow velocity to the lift is more than the influence from the rotor position.The research is of guiding significance to the type selection and design of zero-speed fin stabilizer.

zero-speed；fin Stabilizer；lift characteristics；Transition SST；NACA0012

U664.72

A

1008-830X（2016）02-0150-05

2015-11-23

浙江省自然科學(xué)基金（LY14E090002）；國(guó)家自然科學(xué)基金（51409232）

周衛(wèi)（1987-），男，碩士研究生，研究方向：農(nóng)業(yè)機(jī)械化.