水翼渦激振動的數(shù)值模擬研究

劉胡濤，張懷新，，姚慧嵐

（1.上海交通大學(xué) 船舶海洋與建筑工程學(xué)院，上海 200240；2.高新船舶與深海開發(fā)裝備協(xié)同創(chuàng)新中心（船海協(xié)創(chuàng)中心），上海 200240）

水翼渦激振動的數(shù)值模擬研究

劉胡濤1，張懷新1，2，姚慧嵐1

（1.上海交通大學(xué) 船舶海洋與建筑工程學(xué)院，上海 200240；2.高新船舶與深海開發(fā)裝備協(xié)同創(chuàng)新中心（船海協(xié)創(chuàng)中心），上海 200240）

對二維水翼結(jié)構(gòu)進行流固耦合運動分析，利用大渦模擬方法計算高雷諾數(shù)下水翼繞流場，流場力作用于二自由度剛體上導(dǎo)致周期性的垂蕩和轉(zhuǎn)動，龍格庫塔法求解剛體水翼的運動方程，位移參數(shù)作為下一時間步流場計算的邊界條件，具體通過編譯自定義函數(shù)控制剛體運動和流場網(wǎng)格變化。探討了初始攻角、剛心位置以及來流速度對水翼振動的影響，發(fā)現(xiàn)在來流速度增大至某一數(shù)值時水翼發(fā)生了顫振現(xiàn)象，并在時域與頻域上對顫振的發(fā)生機理進行深入探討。

水翼繞流；流固耦合；顫振

0　引　言

隨著水面艦艇和潛艇航速的提升，水彈性振動對船體強度以及疲勞壽命產(chǎn)生的影響已引起了造船界的廣泛關(guān)注。艦艇附件和海洋結(jié)構(gòu)部件在粘性流體中運動時，引起流動分離，生成旋渦，造成渦激振動，可能導(dǎo)致結(jié)構(gòu)破壞。利用水彈性的研究結(jié)果，分析粘性流體通過翼型振蕩物的振動研究，國內(nèi)外的研究并不多?，F(xiàn)有文獻結(jié)果大多是假定圓柱、機翼靜止不動的情況下進行的。多數(shù)情況下水翼不是完全固定不動的，而是在彈性軸處受到非線性約束，水翼在流體的作用下會產(chǎn)生渦的脫落等復(fù)雜現(xiàn)象，由此流體力作用于結(jié)構(gòu)物本身產(chǎn)生耦合效應(yīng)渦激振動。如水翼艇中的水翼，潛艇中的舵翼等船舶構(gòu)件作為動力學(xué)系統(tǒng)，設(shè)計不當時，在流場中高速運動時將會產(chǎn)生有害的渦激顫振，導(dǎo)致結(jié)構(gòu)破壞。

流固耦合研究早期集中在氣動彈性力學(xué)，發(fā)展至今已經(jīng)有較為成熟的理論體系，其大部分理論結(jié)果都可直接應(yīng)用于水彈性力學(xué)中。早在 20 世紀 30 年代Theodorsen［1］提出基于非定常氣動力線性化精確解的不可壓縮、無粘流的 Theodorsen 非定常氣動力模型，該模型同樣可以處理二維水翼的非定常力。Huang等［2］對NACA0012 翼型進行實驗分析，并探討了翼根連接處與翼梢對渦激振動的影響。Jung［3］通過實驗對低雷諾數(shù)下機翼繞流形成的渦街，發(fā)現(xiàn)機翼振動情況下的渦脫頻率相對于機翼靜止不動時要小。王囡囡［4］采用片條理論的簡化非定常氣動力方法研究機翼氣動彈性系統(tǒng)的動態(tài)特性，提出了一種基于實測柔度的顫振主動控制方法。另一方面至 20 世紀 70 年代后期，船舶水彈性力學(xué)理論和分析方法取得了顯著的進展，并廣泛應(yīng)用到船舶與海洋工程領(lǐng)域的流固耦合分析中。早期的工作主要集中在試驗研究，Davidson 實驗室對水翼做了大量的試驗研究，主要集中在以 2 個自由度理想翼型為出發(fā)點，觀察測量水翼在水洞試驗中的振動、空泡、壓力分布等。Henry［5］進行了一系列的模型試驗，對水翼各參數(shù)以及臨界顫振速度的影響作了理論上探討。余志興［6-7］運用 N-S 方程數(shù)值模擬了粘性流場中二維機翼的水彈性振動，在此基礎(chǔ)上，劉曉宙［8］利用有均流時二維運動物體的聲輻射方程，研究了大攻角下二維水翼的聲輻射，結(jié)果表明：渦脫頻率和機翼的固有振動頻率一致時，聲輻射達到最大。Chae［9-10］對2個自由度下的二維水翼進行了一系列系統(tǒng)的研究，探討了附加水質(zhì)量、水動力阻尼、非線性流固耦合響應(yīng)等對水翼穩(wěn)定性的影響。

本文建立了 2 個自由度下二元水翼的流固耦合模型，分別利用大渦模擬方法和龍格庫塔法求解流體結(jié)構(gòu)方程，界面數(shù)據(jù)傳遞通過自編函數(shù)實現(xiàn)。對NACA0012 翼型在不同工況下的繞流振動進行分析，并分別在時域和頻域上對振動情況進行分析。

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

1.1剛體模型

本文將探討懸臂式、二維翼型在不可壓縮流體中的動態(tài)響應(yīng)與穩(wěn)定邊界，假定水翼僅在翼展方向存在2 個自由度運動，即對應(yīng)船舶六自由度中的垂蕩與縱搖，如圖1 所示。

圖1　兩自由度水翼模型Fig.1　Two degrees of freedom hydrofoil models

其中，水翼弦長 c=2b，來流速度為 U，初始攻角為 α。EA 為水翼彈性中心（也即是二維翼型剛心），位于翼弦中點后 ab 處，a 為彈性中心到翼弦中點距離的無量綱量，b 為半弦長。翼型在彈性軸處由一個線彈簧和一個扭轉(zhuǎn)彈簧支撐在彈性軸處，垂直方向位移與旋轉(zhuǎn)均是在彈性軸上測量，h（t）向上為正，θ（t）逆時針為正，分別稱為水翼的垂蕩運動和縱搖運動。HC 為水動力中心，位于彈性中心前緣 eb 處，即是流體作用力合力的作用點。重心 CG 位于彈性中心后緣 xθb 處，xθ為重心到彈性中心距離的無量綱量。

兩自由度水翼無量綱運動方程為：

上述方程組可簡化為：

式中：Ms，Cs，Ks分別為水翼慣性、阻尼和剛度矩陣；F 為流場力矩陣；X 為兩自由度上廣義位移矩陣。

1.2流體模型

如圖1 所示，流場分為 3 個區(qū)域：水翼近場加密網(wǎng)格為動區(qū)域；遠場采用非結(jié)構(gòu)網(wǎng)格為變形區(qū)域；延長水翼后方網(wǎng)格以觀察尾流。

圖2　水翼計算網(wǎng)格Fig.2　Grid of hydrofoil

基于商業(yè)軟件 Fluent 有限體積法對二維水翼進行瞬態(tài)繞流分析，利用 UDF 控制剛體與網(wǎng)格的運動。大渦模擬相較于平均雷諾數(shù)法可以提供豐富的大渦旋信息，且能計算出壁面處的壓強脈動。本文在對 Fluent中幾種經(jīng)典湍流模型進行計算比較后選用 LES 作為后續(xù)的計算模型。

流場的連續(xù)性方程和不可壓縮 N-S 方程可以寫為：

對于 Smagorinsky-Lilly 模型

1.3流固耦合設(shè)置

建立 CFD 模型和計算結(jié)構(gòu)動力學(xué)模型后，需要建立有效的界面耦合程序來聯(lián)系這 2 個單獨的計算模型。耦合界面不匹配網(wǎng)格間信息的準確交換，直接影響耦合計算精度，是獲得真實物理解的關(guān)鍵環(huán)節(jié)之一。

2　計算結(jié)果及分析

水翼結(jié)構(gòu)參數(shù)取值如下：總弦長 c=1 m，垂蕩、縱搖固有頻率 ωh=ωθ=2s?1，中心對彈性軸的回轉(zhuǎn)半徑 r=0.8，質(zhì)量比 μ=m/（πρb2）=14.0，無量綱系數(shù)a=?1,xθ=0.1。計算中流體與結(jié)構(gòu)的時間步長均取為 0.01 s，是結(jié)構(gòu)自振周期的 0.318%。

對 1 m/s 來流下渦激振動采用大渦模擬數(shù)值分析，剛心位置取在水翼前緣。圖中流場參數(shù)為 Re=106。

圖3　.耦合算法流程圖Fig.3　Coupled algorithm flowsheet

從圖4 可看出，10° 攻角下水翼后緣渦街脫落較為明顯，在水翼前緣和尾翼處均有渦的形成，以平行于來流的方向向尾部脫離，并在尾流中形成渦街。而0°攻角下由于水翼振幅很小，僅在尾翼處有渦的形成，且很快消散。目前關(guān)于翼型繞流的渦激振動研究，為了獲得翼型繞流周期性的渦發(fā)放，大多集中在大攻角工況下。本文以 10°攻角作為主要分析對象，發(fā)現(xiàn)大渦模擬能夠較好的捕捉復(fù)雜湍流中的細節(jié)，而在一般 k-ω等雷諾平均模型中，小攻角水翼繞流無法捕捉到尾部渦這一特征。

保持其他參數(shù)不變，改變初始攻角對比其對水翼繞流振動的影響，垂蕩運動和縱搖運動結(jié)果如圖6 所示。

顯而易見，攻角越大，其升力與力矩就越大，相應(yīng) 2 個自由度上的位移就越大。從圖中可看到，不同攻角下的振蕩平衡位置隨攻角增大而增大。在 0° 攻角下，僅在水翼尾部形成較小的渦，且很快耗散掉，圖5 中可以看出，0° 攻角下，水翼振蕩位移的幅值非常小，且頻率較高，沒有體現(xiàn)出一定的周期性。本文嘗試將 0° 攻角來流速度提升至 4 m/s，2 個自由度上的位移仍保持非常小，且周期性較為紊亂。初步得出結(jié)論：0° 攻角在沒有初始擾動的情況下很難產(chǎn)生大幅振動的情況，攻角越大，振動平衡位置越偏離初始位置。且垂蕩運動幅值衰減較明顯，縱搖運動幅值衰減較慢。

圖4　NACA0012 水翼壓力、流場以及渦量云圖，Re=106Fig.4　Pressure，velocity and vorticity contours of hydrofoil Re=106

圖5　不同攻角下振動情況對比 U=1 m/sFig.5　Hydrofoil vibration under different attack angle

圖6　剛心位置對振動情況的影響 U=1 m/s，α=10°Fig.6　Hydrofoil vibration under different stiffness center

根據(jù)氣動彈性力學(xué)理論，重心位于剛心線之前就不會發(fā)生顫振，飛行器機翼設(shè)計時多采用在翼面前緣防止配重，使重心前移的方式來避免顫振。當重心位置一定時，應(yīng)使剛心靠后一些（與重心靠近），以此來提高顫振速度。對于水翼而言，剛心位置同樣對顫振臨界條件有著至關(guān)重要的作用。通過升力力矩公式可以看出，繞彈性軸力矩的大小與剛心位置有關(guān)，剛心位置的改變直接影響到力矩的大小，進而影響到縱搖運動，由于流體-結(jié)構(gòu)的耦合作用，結(jié)構(gòu)運動狀態(tài)的改變將導(dǎo)致流場變化從而導(dǎo)致垂蕩運動的變化。圖6顯示 2 個自由度上位移都出現(xiàn)了一定幅度的衰減，且剛心位置對縱搖振動平衡位置的影響較垂蕩運動更加明顯。

從圖7可見，隨著來流速度增大，水翼在 2 個自由度上振動不斷加劇，可以看到 1.67 m/s 工況下縱搖自由度振幅隨著時間不斷增大，而垂蕩自由度的振幅經(jīng)過小幅度減小后又繼續(xù)增大直到相對保持平穩(wěn)。初步判定，本文所采用水翼的顫振臨界速度在 1.67 m/s附近。其中，1 m/s 時，升力和力矩在數(shù)值上保持微小的波動，可近似看作是固定力作用在水翼上，由于運動過程中的阻尼作用，振動幅度隨時間逐漸衰減。而 2 m/s時，升力和剛心處力矩幅度隨時間不斷增大，從而導(dǎo)致 2 個自由度上的位移也不斷增大，且此時水翼縱搖自由度瞬時水動力與彈性位移之間有約 180°的位相差，水翼已經(jīng)進入顫振狀態(tài)，振動幅度不斷擴大從而導(dǎo)致結(jié)構(gòu)的破壞。為了更直觀地觀察水翼顫振情況，做出以上 3 組速度下的相圖和對應(yīng)，的譜分析圖（見圖8）。

來流速度較小時，流體對水翼做負功，水動力對水翼運動起阻尼作用，振動逐漸衰減最終穩(wěn)定于某一平衡位置。速度為 1 m/s 時，由相圖分析看出系統(tǒng)未出現(xiàn)極限環(huán)運動，2 個自由度上的振幅不斷減小。由頻譜圖可知，垂蕩和縱搖自由度均由 2 個頻率成分疊加而成，其中垂蕩自由度主要頻率成分為 0.35 Hz、縱搖自由度主要頻率成分為 0.18 Hz 和 0.35 Hz。來流速度為 1.67 m/s 時，此時垂蕩自由度含有 2 個頻率分別為0.195 Hz 和 0.34 Hz?？v搖自由度的主要頻率為 0.195 Hz，由圖8可知，此速度下縱搖方向已不穩(wěn)定，說明頻率較小的成分起到主要作用時系統(tǒng)會發(fā)生不穩(wěn)定現(xiàn)象。速度為 2 m/s 時，2 個自由度上的主要頻率成分均為0.195 Hz，與 1.67 m/s 縱搖不穩(wěn)定時的頻率相同，次要頻率成分為 0.3 Hz 且幅值遠小于 0.195 Hz，此時系統(tǒng)發(fā)生顫振，也同樣印證了上述說法，即較小頻率成分起主要作用時系統(tǒng)不穩(wěn)定。表明來流速度增大至越過臨界速度后，水動力做正功，水翼發(fā)生自激振動，振幅不斷增大而導(dǎo)致振動發(fā)散，顫振是指結(jié)構(gòu)物具有 2個自由度以上的以同一頻率耦合起來的振動現(xiàn)象［11］。由頻譜圖可以看出，隨速度增大，2 個自由度上的主要頻率成分趨于一致，且低頻成分逐漸起到主導(dǎo)作用，另一個頻率成分的影響逐漸降低。

圖7　不同來流速度下振動情況對比，α=10°，a=-1Fig.7　Hydrofoil vibration under different velocity

圖8　不同來流速度下縱搖運動相圖及頻譜Fig.8　Motion phase diagram and velocity spectrum in various velocity

圖9　不同流速下水翼的聲學(xué)指向性圖Fig.9　Acoustic directivity chart in various velocity

以水翼壁面上的脈動壓力作為偶極子噪聲源，對水翼渦激振動造成的水動力噪聲進行分析，求得圓心為翼弦中點，半徑為 1.5 m 圓上特征點的聲壓值，作出水翼流噪聲的聲學(xué)指向性圖。從圖9 可見，由 1 m/s增至 2 m/s 流噪聲增幅較大，而后增速減緩，這是由于顫振導(dǎo)致振幅的不斷增大導(dǎo)致的。且流噪聲在水翼首、尾位置處較小，這是因為來流和尾流對噪聲起著“掩蔽效應(yīng)”，來流和尾流在靠近壁面處壓力變化梯度大，流場內(nèi)的不均勻波動會抵消一部分流噪聲的傳播，因此在這一區(qū)域流噪聲數(shù)值相對較小。圖10 為特征點（2.0，0）處水翼低頻段聲壓級對比圖?？梢钥闯?，流噪聲主要集中在低頻區(qū)域，隨著速度增大，高頻成分逐漸占據(jù)主要作用，如 2 m/s 時，聲級主要集中在 0～5 Hz 之間，在 3 Hz 處出現(xiàn)峰值，速度增至 4 m/s時，聲級則集中在 0～20 Hz 之間，并出現(xiàn)了多個峰值。聯(lián)合頻譜分析可看出，振動頻率與噪聲頻率都隨來流速度增大而增大。

圖10　不同流速下特征點聲壓值對比圖Fig.10　SPL of feature points in various velocity

3　結(jié)　論

本文通過數(shù)值方法分析二元水翼在紊流中的振動情況，采用弱耦合的方式聯(lián)立大渦模擬湍流模型和剛體運動方程求解了水翼在 2 個自由度上的運動，在此基礎(chǔ)上分別對比了初始攻角、剛心位置以及來流速度對振動的影響，得出以下幾點結(jié)論：

1）大渦模擬方法能夠較好的捕捉高雷諾數(shù)、復(fù)雜湍流中的細節(jié)，模擬出小攻角水翼繞流后完整的渦街現(xiàn)象。

2）水翼振動狀態(tài)對系統(tǒng)的初始值具有依賴性，在沒有擾動的情況下，水翼在 0° 攻角時始終保持微幅振動，攻角越大，振動平衡位置越偏離初始位置，且垂蕩運動幅值衰減較明顯。

3）航行速度對水翼振動影響較為直觀，且隨著速度的不斷增大，水翼出現(xiàn)顫振現(xiàn)象，即振幅不衰減的自激振動。雖然水中顫振產(chǎn)生的條件比較苛刻，但是隨流速增大而產(chǎn)生的振動加劇現(xiàn)象仍不容忽視。

4）流噪聲隨來流速度增大而增大，其在低頻段的數(shù)值高于高頻段的數(shù)值，它能夠直觀地反應(yīng)出結(jié)構(gòu)振動對流場的擾動情況。

［1］THEODORSEN T，MUTCHLER W H.General theory of aerodynamic instability and the mechanism of flutter［R］.NACA Report No.496，NACA，1935.

［2］HUANG R F，LIN C L.Vortex shedding and shear-layer instability of wing at low-Reynolds numbers［J］.AIAA Journal，1995，33（8）: 1398-1403.

［3］JUNG Y W，PARK S O.Vortex-shedding characteristics in the wake of an oscillating airfoil at low Reynolds number［J］.Journal of Fluids and Structures，2005，20（3）: 451-464.

［4］王囡囡.二元機翼顫振及其主動控制的研究［D］.徐州: 中國礦業(yè)大學(xué)，2013.WANG Nan-nan.Research on flutter and active control for two-dimensional airfoil［D］.Xuzhou: China University of Mining and Technology，2013.

［5］HENRY C J.Hydrofoil flutter phenomenon and airfoil flutter theory［M］.New Jersey: Davidson Laboratory，1961.

［6］余志興.粘性流場中的水彈性計算［D］.上海: 上海交通大學(xué)，1999.YU Zhi-xing.Hydroelastic dynamics in viscous flow［D］.Shanghai: Shanghai Jiao Tong University，1999.

［7］余志興，劉應(yīng)中，繆國平.二維機翼彈簧系統(tǒng)的渦激振動［J］.船舶力學(xué)，2002，6（5）: 25-32.YU Zhi-xing，LIU Ying-zhong，MIAO Guo-ping.Vortex-induced vibration of two-dimensional wing-spring coupled system［J］.Journal of Ship Mechanics，2002，6（5）: 25-32.

［8］劉曉宙，繆國平，余志新，等.流體通過渦激振動機翼的聲輻射研究［J］.聲學(xué)學(xué)報，2005，30（1）: 55-62.LIU Xiao-zhou，MIAO Guo-ping，YU Zhi-xin，et al.Research on the sound generation by flow around vortex-induced vibrated aerofoil［J］.Acta Acustica，2005，30（1）: 55-62.

［9］CHAE E J，AKCABAY D T，YOUNG Y L.Dynamic response and stability of a flapping foil in a dense and viscous fluid［J］.Physics of Fluids，2013，25（10）: 104106.

［10］CHAE E J.Dynamic response and stability of flexible hydrofoils in incompressible and viscous flow［D］.Michigan: University of Michigan，2015.

［11］程貫一，王寶壽，張效慈.水彈性力學(xué)——基本原理與工程應(yīng)用［M］上海: 上海交通大學(xué)出版社，2013.CHENG Guan-yi，WANG Bao-shou，ZHANG Xiao-ci.Hydroelasticity——the basics with applications［M］.Shanghai: Shanghai Jiao Tong University，2013.

Research on vortex induced vibration of hydrofoil

LIU Hu-tao1，ZHANG Huai-xin1，2，YAO Hui-lan1
（1.School of Naval Architecture，Ocean and Civil Engineering，Shanghai Jiaotong University，Shanghai 200240，China；2.Collaborative Innovation Center for Advanced Ship and Deep-Sea Exploration（CISSE），Shanghai 200240，China）

The objective of this research is to investigate the fluid-structure interaction motion of a two-dimensional hydrofoil in viscous flow.Large-eddy simulation method is used to calculate the flow field around hydrofoil，and Runge-Kutta method is applied to solve equations of rigid body motion.Flow results and motion results exchange data at the hydrofoil wall surface，by compiling a user-defined functions to control rigid body motion and flow field grid update.The influence of the hydrofoil parameters，i.e.attack angle，elastic axis position and flow velocity，on vibration are discussed.Once the velocity increase to a certain value，hydrofoil flutter phenomenon occurred.Also，discussion on time and frequency domain is held to investigate the flutter phenomenon.

two-dimensional wing flow；fluid Structure Interaction；flutter

U661.1

A

1672-7619（2016）06-0007-07

10.3404/j.issn.1672-7619.2016.06.002

2015-12-09；

2016-03-17

國家自然科學(xué)基金資助項目 （11272213，51479116）

劉胡濤（1991-），男，碩士研究生，研究方向為船舶水動力學(xué)。