低質(zhì)量比舵系統(tǒng)流致振動(dòng)特征研究

王人鳳， 尤云祥， 陳 科， 胡曉峰

(上海交通大學(xué) 海洋工程國(guó)家重點(diǎn)實(shí)驗(yàn)室 高新船舶與深海開發(fā)裝備協(xié)同創(chuàng)新中心，上海 200240)

顫振是一種由彈性力、慣性力、阻尼力和流體自激力共同作用引起的彈性不穩(wěn)定現(xiàn)象。由于這類振動(dòng)有可能演化為發(fā)散振動(dòng)，故該現(xiàn)象一旦發(fā)生，結(jié)構(gòu)就存在損壞的風(fēng)險(xiǎn)。經(jīng)典顫振理論認(rèn)為質(zhì)量比(μ=m/πρb2l，m為系統(tǒng)質(zhì)量，ρ為流體密度，b為舵的半弦長(zhǎng)，l為舵的展長(zhǎng))低于臨界質(zhì)量比時(shí)顫振便不會(huì)發(fā)生，但一些試驗(yàn)表明，當(dāng)水下升力體的質(zhì)量比較小時(shí)也有可能發(fā)生顫振[1]。顫振理論最初應(yīng)用于航空領(lǐng)域，后又被推廣到船舶領(lǐng)域。機(jī)翼的臨界顫振速度隨著機(jī)翼質(zhì)量、質(zhì)心到剛心的距離以及結(jié)構(gòu)阻尼系數(shù)的增加而增大。當(dāng)結(jié)構(gòu)具有立方非線性剛度時(shí)，系統(tǒng)會(huì)在一定的速度范圍內(nèi)做穩(wěn)定的極限環(huán)運(yùn)動(dòng)，并且極限環(huán)的幅值會(huì)隨著流速的增加而增大。但是，當(dāng)流速增大到一定程度時(shí)，系統(tǒng)進(jìn)行發(fā)散振動(dòng)。在求解振動(dòng)方程時(shí)，不同的初始積分條件也會(huì)對(duì)系統(tǒng)極限環(huán)振動(dòng)的幅值產(chǎn)生影響[2]。

張琪昌等[3]針對(duì)含立方非線性剛度二元機(jī)翼系統(tǒng)顫振時(shí)的極限環(huán)振動(dòng)開展了相關(guān)研究，對(duì)分岔點(diǎn)類別和平衡點(diǎn)的性質(zhì)進(jìn)行了定性的分析，進(jìn)而研究了系統(tǒng)參數(shù)尤其是機(jī)翼外形變化對(duì)極限環(huán)顫振穩(wěn)定性和幅值的影響。發(fā)現(xiàn)當(dāng)平衡點(diǎn)為穩(wěn)定的焦點(diǎn)時(shí)，機(jī)翼會(huì)進(jìn)行衰減振動(dòng)，滿足設(shè)計(jì)要求。但是，在平衡點(diǎn)的附近也有可能出現(xiàn)穩(wěn)定的極限環(huán)顫振運(yùn)動(dòng)，甚至是不穩(wěn)定的極限環(huán)顫振運(yùn)動(dòng)。在工程上，對(duì)于小振幅穩(wěn)定極限環(huán)顫振，只要不引起疲勞破壞就可以滿足要求；而對(duì)于不穩(wěn)定的極限環(huán)顫振運(yùn)動(dòng)，必須避免。設(shè)計(jì)中可以通過改變系統(tǒng)的參數(shù)來抑制振幅和臨界顫振速度的大小，從而避免不穩(wěn)定極限環(huán)顫振運(yùn)動(dòng)的產(chǎn)生。趙海等[4-5]求解某彈性系統(tǒng)的顫振速度并進(jìn)行了數(shù)值驗(yàn)證，分析了系統(tǒng)各參數(shù)對(duì)顫振速度的影響，并討論了非線性剛度系數(shù)對(duì)系統(tǒng)沉浮和俯仰自由度顫振極限環(huán)幅值的影響。結(jié)果表明，忽略非線性阻尼項(xiàng)會(huì)對(duì)超聲速流二元機(jī)翼的顫振計(jì)算造成誤差。其研究還顯示，雖然系統(tǒng)運(yùn)動(dòng)微分方程具有對(duì)稱性，但兩個(gè)自由度上的結(jié)構(gòu)參數(shù)對(duì)系統(tǒng)臨界顫振速度的影響卻不具有對(duì)稱性。張瑜等[6-7]利用變步長(zhǎng)四階龍格庫(kù)塔法對(duì)二元機(jī)翼系統(tǒng)的運(yùn)動(dòng)微分方程進(jìn)行數(shù)值模擬，對(duì)不同流速區(qū)間內(nèi)平衡點(diǎn)、極限環(huán)個(gè)數(shù)以及穩(wěn)定性進(jìn)行研究。同時(shí)還對(duì)極限環(huán)分叉問題進(jìn)行了數(shù)值模擬，結(jié)果表明系統(tǒng)的極限環(huán)也有可能產(chǎn)生和平衡點(diǎn)一樣的分岔現(xiàn)象，且隨著無量綱流速的增大，極限環(huán)振動(dòng)會(huì)經(jīng)過倍周期分叉轉(zhuǎn)化為混沌運(yùn)動(dòng)。

Hamdani等[8-9]將有限元方法應(yīng)用于翼型的動(dòng)力特性分析中，一方面可以得到翼型表面的壓力分布，另一方面也可以得到翼型附近的流體分布。翼型在低雷諾數(shù)下同樣會(huì)出現(xiàn)運(yùn)動(dòng)速度突然增大或減小的情況，并且在此期間會(huì)產(chǎn)生一個(gè)較大的升力，導(dǎo)致翼形的進(jìn)行明顯的沉浮和俯仰運(yùn)動(dòng)。在這種不穩(wěn)定的運(yùn)動(dòng)發(fā)生期間，會(huì)在翼型的上下表面不斷出現(xiàn)一層層范圍更大的漩渦。Münch等[10]利用CFX模擬了NACA0009翼型在湍流中的動(dòng)力特性。翼型的初始攻角為2°，流速為5～15 m/s，計(jì)算結(jié)果與實(shí)驗(yàn)數(shù)據(jù)吻合，誤差小于2%。Gnesin等[11]利用一種考慮旋轉(zhuǎn)葉片振動(dòng)因素的三維流固耦合求解方法研究了渦輪轉(zhuǎn)子的振動(dòng)問題。發(fā)現(xiàn)在葉片振動(dòng)的幅值-頻率譜上存在高頻和低頻，但是這些頻率并不是轉(zhuǎn)子的倍頻。Amiralaei等[12]為分析NACA0012翼型在低雷諾數(shù)流體中的振動(dòng)特性，分別選取攻角2°～10°，使用有限體積法求解Navier-Stokes方程，并將計(jì)算結(jié)果中的瞬時(shí)升力系數(shù)與利用Theodorsen理論求解的結(jié)果進(jìn)行比較，發(fā)現(xiàn)雷諾數(shù)不僅對(duì)翼型的升阻力系數(shù)有影響，對(duì)翼型的振動(dòng)特征也存在明顯的影響。Akcabay等[13]利用CFX求解非定常雷諾平均Navier-Stokes方程分析空泡對(duì)NACA66水翼水彈性穩(wěn)定性的影響。隨著空泡量的增加，翼型的平均升力隨之減小而平均阻力增大，同時(shí)其平均變形幅度也隨之減小。此外，由流動(dòng)引起的振動(dòng)幅值和阻力隨著翼型剛度的減小而增大。劉胡濤等[14]利用流固耦合方法對(duì)二元水翼的運(yùn)動(dòng)進(jìn)行數(shù)值模擬，分析了初始攻角、剛心位置以及來流速度對(duì)水翼振動(dòng)特性的影響。使用大渦模擬方法計(jì)算高雷諾數(shù)下水翼繞流場(chǎng)以及流場(chǎng)力作用于兩自由度剛體上所導(dǎo)致的周期性沉浮和俯仰運(yùn)動(dòng)。采用龍格庫(kù)塔法求解剛體水翼的運(yùn)動(dòng)方程，位移參數(shù)作為下一時(shí)間步流場(chǎng)計(jì)算的邊界條件。結(jié)果表明，水翼振動(dòng)狀態(tài)對(duì)系統(tǒng)的初始值具有依賴性，在沒有擾動(dòng)的情況下，水翼在零攻角時(shí)始終保持微幅振動(dòng)。隨著初始攻角的增大，振動(dòng)平衡位置越偏離初始位置，且沉浮運(yùn)動(dòng)幅值衰減較明顯。來流速度對(duì)水翼振動(dòng)影響較為顯著，并且當(dāng)來流速度增大到一定程度時(shí)，水翼會(huì)出現(xiàn)顫振現(xiàn)象。此外，還發(fā)現(xiàn)水翼顫振產(chǎn)生的條件比較苛刻。盡管如此，隨著流速的增大，系統(tǒng)產(chǎn)生的振動(dòng)加劇現(xiàn)象仍不容忽視。

若舵具有較小展弦比，利用兩自由度系統(tǒng)動(dòng)力方程分析振動(dòng)特性時(shí)需要考慮其三維效應(yīng)，所以有必要對(duì)理論模型進(jìn)行相應(yīng)的參數(shù)修正。此外，鑒于低質(zhì)量比(μ<1)舵系統(tǒng)顫振相關(guān)研究較少，通過數(shù)值模擬方法分析顫振發(fā)生時(shí)舵周圍的流場(chǎng)變化也很有必要。

1 模型和計(jì)算方法

1.1 小展弦比兩自由度系統(tǒng)動(dòng)力模型

舵系統(tǒng)由一個(gè)截面為NACA0017的舵、一根剛性軸以及與其連接的彈性裝置組成，剛性軸與舵固接。其中，舵的展長(zhǎng)l=0.39 m，兩端弦長(zhǎng)不同，舵根側(cè)面(較大的側(cè)面)的弦長(zhǎng)為0.3 m，舵稍側(cè)面(較小的側(cè)面)的弦長(zhǎng)為0.2 m，則中截面的弦長(zhǎng)c=0.25 m?？梢詫⒍嫦到y(tǒng)簡(jiǎn)化為經(jīng)典兩自由度舵模型，如圖1所示。其主要結(jié)構(gòu)參數(shù)如表1所示。

圖1 有限展長(zhǎng)舵的兩自由度振動(dòng)理論計(jì)算模型Fig.1 Two-degree-freedom vibration mechanical model for a foil with a finite span

m/kgμAOA/(°)kh/(N·m-1)kα/(N·m·rad-1)xα/(N·m-1)a7.60.451.4×1062820.04-0.38

當(dāng)舵做簡(jiǎn)諧運(yùn)動(dòng)時(shí)，會(huì)受到升力L(向上為正)和俯仰力矩Tα(迎流抬頭為正)的共同作用，并在流體動(dòng)力的激勵(lì)下發(fā)生兩個(gè)自由度的運(yùn)動(dòng)，一個(gè)是舵-軸部件的沉浮運(yùn)動(dòng)，位移為h，向下為正；另一個(gè)是舵隨剛性軸轉(zhuǎn)動(dòng)而產(chǎn)生的俯仰運(yùn)動(dòng)，俯仰角度為α，迎流抬頭為正。

(1)

式中：V∞為流速；C(k)為Theodorsen函數(shù)； 式中與C(k)無關(guān)的項(xiàng)來自慣性效應(yīng)。

考慮展弦比對(duì)結(jié)構(gòu)振動(dòng)的影響，對(duì)該式進(jìn)行修正。代表慣性效應(yīng)的項(xiàng)添加附加質(zhì)量修正系數(shù)ε，而環(huán)量項(xiàng)添加環(huán)量修正系數(shù)δ，得到

(2)

其中，

(3)

式中：AR=l/(ccosAOA)，AOA為初始攻角；τ為舵的形狀參數(shù)。

通過Fourier變換和引入Wagner函數(shù)[15]，結(jié)合舵的兩自由度任意運(yùn)動(dòng)微分方程

(4)

可以得到

(5)

其中，

(6)

將式(5)寫為矩陣形式，有

(7)

(8)

從而得到兩自由度系統(tǒng)任意運(yùn)動(dòng)的無量綱動(dòng)力方程為

(9)

1.2 流體模型

使用CFD(Computational fluid dynamics)商業(yè)軟件中的彈簧構(gòu)件生成系統(tǒng)特定的支撐剛度，并利用UDF(User Defined Function)通過添加與俯仰運(yùn)動(dòng)方向相反的俯仰力矩來模擬相應(yīng)的扭轉(zhuǎn)剛度。對(duì)于流體域，左側(cè)以及前后壁面為速度入口，上下為固壁面，右側(cè)為壓力出口，如圖2所示。

圖2 數(shù)值模擬幾何圖Fig.2 Geometric model of numerical simulation

當(dāng)?shù)退俸叫袝r(shí)，舵系統(tǒng)周圍的流體可以視作不可壓流體。質(zhì)量守恒方程和動(dòng)量守恒方程組成了黏性不可壓流體的基本控制方程，而控制方程正是這些守恒定律的數(shù)學(xué)描述?，F(xiàn)實(shí)中的湍流流動(dòng)是一個(gè)隨機(jī)變化的脈動(dòng)量，由于很難在計(jì)算中完全反映這些脈動(dòng)量，所以引入雷諾時(shí)均的概念，將任一流動(dòng)變量的瞬時(shí)值進(jìn)行分解，并看作平均值和脈動(dòng)值之和。

三維瞬態(tài)不可壓流體質(zhì)量守恒方程為

(10)

時(shí)均值的動(dòng)量守恒方程，即雷諾時(shí)均納維爾—斯托克斯(Reynolds Average Navier-Stokes, RANS)為[16]

(11)

式中： (u,v,w)為流速在直角坐標(biāo)系中的三個(gè)分量;ρ為流體的密度;p為壓力;μ為流體的動(dòng)力黏性系數(shù)；g為重力加速度。

由Boussinesq假設(shè)，引入湍流動(dòng)力黏性系數(shù)μt，可以得到雷諾應(yīng)力張量和流體時(shí)均速度梯度之間的關(guān)系

只要選擇合適的湍流模型將μt確定，就可以使RANS封閉。本文選擇的SST(Shear Stress Transport)k-ω模型(經(jīng)剪切應(yīng)力修正的兩方程模型)結(jié)合了標(biāo)準(zhǔn)兩方程模型k-ω和k-ε的優(yōu)點(diǎn)，在近壁面區(qū)域利用k-ω模型，而在遠(yuǎn)離附面層的區(qū)域使用k-ε模型，可以精確地計(jì)算物體表面流體的逆壓梯度以得到物體附近流體的細(xì)微變化[17-18]。

1.3 流固耦合設(shè)置

把舵系統(tǒng)的運(yùn)動(dòng)看作是一種剛體運(yùn)動(dòng)，其重心位于剛體運(yùn)動(dòng)局部坐標(biāo)系的原點(diǎn)。為了保持和實(shí)際情況一致，在剛心的垂向位置添加了彈性系數(shù)為kh的彈簧組構(gòu)件，并通過UDF設(shè)置了舵系統(tǒng)繞剛性軸方向的俯仰力矩，其表達(dá)式為

Tα=-kαθ

(12)

式中：θ為舵轉(zhuǎn)動(dòng)的弧度。

舵系統(tǒng)在運(yùn)動(dòng)過程中所受到的回復(fù)力F和俯仰力矩Tα為

(13)

(14)

式中： [τ]和p為作用的剪應(yīng)力和壓力；r-rc為舵表面任意網(wǎng)格中心到舵重心的距離(用位置矢量表示)；n為舵表面的外法線方向。

得到舵系統(tǒng)的受力之后，需要進(jìn)一步計(jì)算其運(yùn)動(dòng)姿態(tài)。引入舵系統(tǒng)的兩自由度運(yùn)動(dòng)方程可得

(15)

(16)

式中：v為舵系統(tǒng)沉浮運(yùn)動(dòng)的速度；M為慣性矩張量在繞軸方向上的分量；ω為舵系統(tǒng)俯仰運(yùn)動(dòng)角速度。由式(15)和式(16)可以求得舵系統(tǒng)運(yùn)動(dòng)的速度、角速度以及位移，進(jìn)而確定其運(yùn)動(dòng)姿態(tài)。在求解過程中，每個(gè)時(shí)間步長(zhǎng)內(nèi)對(duì)離散方程進(jìn)行多次迭代，當(dāng)滿足收斂條件時(shí)輸出速度和角速度等參數(shù)，并即時(shí)更新網(wǎng)格節(jié)點(diǎn)的位置，從而實(shí)現(xiàn)舵系統(tǒng)在流體中運(yùn)動(dòng)的數(shù)值模擬。

1.4 網(wǎng)格劃分

在舵和流體耦合的過程中，網(wǎng)格移動(dòng)是耦合界面數(shù)據(jù)傳遞的關(guān)鍵。當(dāng)結(jié)構(gòu)的運(yùn)動(dòng)形態(tài)存在較大變化時(shí)，網(wǎng)格移動(dòng)速度和準(zhǔn)確度會(huì)對(duì)數(shù)值模擬結(jié)果的效率和精度產(chǎn)生巨大影響。重疊網(wǎng)格技術(shù)(Overlapping Grid Method，OGM)的使結(jié)構(gòu)網(wǎng)格的自動(dòng)更新成為可能。OGM將覆蓋整個(gè)計(jì)算域的母網(wǎng)格和運(yùn)動(dòng)物體的子網(wǎng)格進(jìn)行疊加來描述物體間的相對(duì)運(yùn)動(dòng)，各區(qū)域中的計(jì)算網(wǎng)格獨(dú)立生成，流場(chǎng)信息通過插值函數(shù)在重疊區(qū)邊界進(jìn)行匹配和耦合，從而實(shí)現(xiàn)網(wǎng)格的自動(dòng)更新。

依據(jù)舵模型的尺寸建立相應(yīng)的幾何模型，并劃分面網(wǎng)格，網(wǎng)格數(shù)為23 834個(gè)，如圖3所示。

圖3 舵的概念、實(shí)物、幾何建模和網(wǎng)格圖Fig.3 Conceptual graph, photograph, geometric model and grid of the hydrofoil

圖4為體網(wǎng)格的劃分，網(wǎng)格數(shù)為122萬個(gè)。其中，對(duì)舵表面附近的網(wǎng)格進(jìn)行了加密。在計(jì)算域內(nèi)采用動(dòng)網(wǎng)格技術(shù)，舵附近和外邊界分別設(shè)置了不同的網(wǎng)格密度，舵附近的網(wǎng)格最密，而越接近外邊界網(wǎng)格越稀疏。重疊網(wǎng)格實(shí)現(xiàn)了不同網(wǎng)格區(qū)域的數(shù)據(jù)傳遞，而通過設(shè)置疏密程度不同的網(wǎng)格，既可以精確反映舵表面附近的流場(chǎng)變化，又可以提高計(jì)算速度。

圖4 網(wǎng)格劃分Fig.4 Mesh generation

2 結(jié)果與分析

2.1 振動(dòng)特性分析

圖5為舵系統(tǒng)振動(dòng)理論值的V-g曲線，在給定的系統(tǒng)參數(shù)下存在一個(gè)臨界航速VF，當(dāng)V∞VF時(shí)，沉浮分支所對(duì)應(yīng)的g為負(fù)值，這時(shí)系統(tǒng)進(jìn)行發(fā)散振動(dòng)。

圖5 舵系統(tǒng)振動(dòng)計(jì)算值V-g曲線圖Fig.5 Computed V-g curve for the vibration of the hydrofoil system

舵系統(tǒng)沉浮和俯仰振動(dòng)幅度的實(shí)驗(yàn)值變化曲線，如圖6所示。系統(tǒng)的兩種振動(dòng)成份均隨著V∞增加而增大。沉浮振動(dòng)幅度曲線在V∞=2.37 m/s時(shí)出現(xiàn)一個(gè)明顯的波動(dòng)，并且此流速所對(duì)應(yīng)的頻譜圖上在原有頻率成份的波峰1附近出現(xiàn)了明顯的伴隨峰值?？梢哉J(rèn)定此時(shí)系統(tǒng)發(fā)生顫振，但是這種顫振與高質(zhì)量比(μ?1)系統(tǒng)的顫振不同，其振動(dòng)V-g曲線向下凸。從實(shí)驗(yàn)現(xiàn)象來看，一方面由于阻尼的存在，另一方面V∞較低并且系統(tǒng)的質(zhì)量較小，所以并沒有出現(xiàn)高速水翼以及機(jī)翼經(jīng)歷顫振時(shí)的結(jié)構(gòu)損毀。

圖6 舵系統(tǒng)振動(dòng)幅度變化和振動(dòng)頻譜圖Fig.6 Amplitude growth of the vibrations and frequency spectrums of the hydrofoil system

圖7為V∞=1.0 m/s時(shí)舵系統(tǒng)振動(dòng)時(shí)歷曲線(見圖7(a)和圖7(c))和相軌線(見圖7(b)和圖7(d))?？梢?，系統(tǒng)的運(yùn)動(dòng)中存在著沉浮和俯仰兩個(gè)成份，分別是系統(tǒng)對(duì)支撐彈簧回復(fù)力和扭轉(zhuǎn)彈簧回復(fù)力矩的響應(yīng)。由于V∞=1.0 m/s沒有達(dá)到系統(tǒng)的臨界顫振速度VF，所以系統(tǒng)進(jìn)行衰減振動(dòng)。系統(tǒng)的動(dòng)能會(huì)在運(yùn)動(dòng)的過程中發(fā)生耗散，表現(xiàn)為系統(tǒng)沉浮和俯仰運(yùn)動(dòng)的幅度逐漸變小，直至達(dá)到穩(wěn)定狀態(tài)。從系統(tǒng)運(yùn)動(dòng)相軌線可以看出，沉浮成份表現(xiàn)為圍繞原點(diǎn)做周期運(yùn)動(dòng)，由于振幅單調(diào)減小，所以相軌線不斷持續(xù)向原點(diǎn)靠近，最終在原點(diǎn)達(dá)到穩(wěn)定狀態(tài)。而俯仰成份的周期運(yùn)動(dòng)圍繞著一個(gè)飄移的中心進(jìn)行。從計(jì)算結(jié)果來看，舵系統(tǒng)流致運(yùn)動(dòng)的沉浮成份是一種簡(jiǎn)諧振動(dòng)，而俯仰成份為兩個(gè)頻率的耦合振動(dòng)。

圖7 V∞=1.0 m/s時(shí)，系統(tǒng)振動(dòng)時(shí)歷曲線和相軌線Fig.7 Time history cures and phase orbits of the hydrofoil system at V∞=1.0 m/s

經(jīng)計(jì)算，舵系統(tǒng)的臨界顫振速度為VF=2.37 m/s，圖8(a)和圖8(c)為系統(tǒng)發(fā)生顫振時(shí)的時(shí)歷曲線圖?？梢姡粮〕煞葑龅确?jiǎn)諧振動(dòng)，即運(yùn)動(dòng)的過程中不會(huì)發(fā)生能量耗散；而俯仰成份在經(jīng)歷了短暫的飄移振動(dòng)之后，也進(jìn)行和沉浮成份相似的等幅振動(dòng)。這是顫振理論中典型的振動(dòng)形態(tài)，在忽略阻尼的作用時(shí)，系統(tǒng)將進(jìn)行無消減的等幅振動(dòng)，如果航行體在高速運(yùn)動(dòng)中發(fā)生這種振動(dòng)有可能造成結(jié)構(gòu)的扭曲、損壞，還有可能產(chǎn)生較大的噪聲導(dǎo)致降低航行體的性能。系統(tǒng)此時(shí)相應(yīng)的運(yùn)動(dòng)相軌線，如圖8(b)和圖8(d)所示。沉浮振動(dòng)和俯仰振動(dòng)最終都演化為典型的極限環(huán)振動(dòng)。因?yàn)橄到y(tǒng)在運(yùn)動(dòng)中不存在能量的耗散，所以能量在動(dòng)能和勢(shì)能之間不斷轉(zhuǎn)化。

圖8 V∞=2.37 m/s時(shí)，系統(tǒng)振動(dòng)時(shí)歷曲線和相軌線Fig.8 Time history cures and phase orbits of the hydrofoil system at V∞=2.37 m/s

圖9為利用未經(jīng)參數(shù)修正的運(yùn)動(dòng)方程求解所得的非等截面系統(tǒng)發(fā)生顫振時(shí)的時(shí)歷曲線圖和相軌線。可見，系統(tǒng)運(yùn)動(dòng)的沉浮成份和俯仰成份都進(jìn)行發(fā)散振動(dòng)。沉浮成份以原點(diǎn)為中心做周期運(yùn)動(dòng)，但是與低流速時(shí)不同，運(yùn)動(dòng)軌跡不斷向外擴(kuò)散；而俯仰成份的運(yùn)動(dòng)軌跡在出現(xiàn)短暫的飄移后，也轉(zhuǎn)化以原點(diǎn)為中心向外不斷做周期性的擴(kuò)散。理論上，出現(xiàn)這種情況是由于系統(tǒng)在運(yùn)動(dòng)的過程中不但沒有發(fā)生能量耗散，還持續(xù)不斷地從流體中吸收能量。這些能量由系統(tǒng)的動(dòng)能轉(zhuǎn)化為勢(shì)能，又由勢(shì)能轉(zhuǎn)化為動(dòng)能，如此循環(huán)往復(fù)便出現(xiàn)了系統(tǒng)運(yùn)動(dòng)過程中沉浮和俯仰成份振幅不斷擴(kuò)大的現(xiàn)象。與圖8進(jìn)行比較，可知在未進(jìn)行參數(shù)修正的前提下，利用兩自由度系統(tǒng)任意運(yùn)動(dòng)方程求解會(huì)導(dǎo)致VF偏小。

當(dāng)來流速度V∞=5.0 m/s時(shí)，非等截面舵系統(tǒng)振動(dòng)的時(shí)歷曲線和相軌線，如圖10所示?？梢?，系統(tǒng)振動(dòng)的沉浮成份和俯仰成份都呈現(xiàn)發(fā)散狀態(tài)，并且幅值迅速擴(kuò)大。圖10(b)和圖10(d)為系統(tǒng)運(yùn)動(dòng)的相軌線圖，沉浮成份和俯仰成份都以原點(diǎn)為中心做周期運(yùn)動(dòng)，運(yùn)動(dòng)軌跡不斷向外擴(kuò)散，并且擴(kuò)散的速度不斷增大。如果不考慮阻尼的影響，現(xiàn)實(shí)中的舵系統(tǒng)若在此流速下進(jìn)行發(fā)散振動(dòng)，則很可能因?yàn)槌粮『透┭稣穹募ぴ龆獾綋p壞。

圖9 V∞=2.37 m/s時(shí)，利用未修正方程求解的時(shí)歷曲線和相軌線Fig.9 Time history cures and phase orbits of the hydrofoil system via unmodified equation for the vibrations at V∞=2.37 m/s

圖10 V∞=5.0 m/s時(shí)，系統(tǒng)振動(dòng)時(shí)歷曲線和相軌線Fig.10 Time history cures and phase orbits of the hydrofoil system at V∞=5.0 m/s

2.2 繞流場(chǎng)分析

圖11為V∞=1.0 m/s時(shí)舵附近的流場(chǎng)變化圖。在運(yùn)動(dòng)最初(見圖11(a))，舵的前緣出現(xiàn)零星的細(xì)小渦，但是并不明顯，舵的后緣附近出現(xiàn)了細(xì)小的尾渦。t=0.05 s(見圖11(b))時(shí)，在舵的前緣出現(xiàn)了明顯的渦量變化，同時(shí)舵的后緣出現(xiàn)明顯脫落渦。t=0.1 s(見圖11(c))時(shí)，在舵上表面從前緣到剛心(約1/4弦長(zhǎng)處)之間出現(xiàn)較大的渦量，而后3/4的上表面有形成多個(gè)細(xì)小的渦。t=0.15 s時(shí)(見圖11(d))，在舵的上表面從剛心到后緣出現(xiàn)了明顯的渦層。t=0.5 s(見圖11(e))時(shí)，在舵的上表面，前緣和剛心之間的渦量繼續(xù)增大，剛心和后緣之間的渦層不斷擴(kuò)大；在舵的下表面，也有出現(xiàn)渦層的趨勢(shì)。由于系統(tǒng)的運(yùn)動(dòng)發(fā)生衰減，所以t=4 s(見圖11(f))時(shí)舵趨于穩(wěn)定，俯仰運(yùn)動(dòng)幅度銳減導(dǎo)致其尾部脫落渦消失，而沉浮運(yùn)動(dòng)幅度銳減導(dǎo)致渦集中在舵的后半部附近。

圖11 V∞=1.0 m/s時(shí)，舵附近流場(chǎng)以及渦量云圖Fig.11 Vorticity contours of the hydrofoil at V∞=1.0 m/s

圖12為V∞=VF=2.37 m/s時(shí)舵附近的流場(chǎng)變化圖?？梢?，在運(yùn)動(dòng)最初(見圖12(a))，和衰減運(yùn)動(dòng)時(shí)相似，舵的后緣附近出現(xiàn)細(xì)小的尾渦，同時(shí)前緣也出現(xiàn)了零星的渦。t=0.05 s(見圖12(b))時(shí)，舵的上表面前緣和半弦長(zhǎng)之間出現(xiàn)了明顯的渦量變化，并且剛心位置的渦最為突出，同時(shí)在尾部出現(xiàn)一系列明顯的脫落渦。t=0.1 s(見圖12(c))時(shí)，在舵的上表面，前緣和剛心之間的渦不斷向后緣擴(kuò)散，并在剛心和后緣之間形成了明顯的渦層。此時(shí)舵的上表面形成一個(gè)完整的渦層，這與衰減運(yùn)動(dòng)時(shí)不同；在舵的下表面，在剛心和后緣之間也有形成渦層的趨勢(shì)。t=0.15 s(見圖12(d))時(shí)，舵上表面的渦層繼續(xù)向后擴(kuò)散，并與尾渦融合。t=2 s(見圖12(e))時(shí)，舵出現(xiàn)大角度的俯仰運(yùn)動(dòng)，此時(shí)前緣和剛心之間出現(xiàn)加大的渦量，而剛心和后緣之間的渦層已經(jīng)完全消失。t=4 s(見圖12(f))時(shí)，只存在前緣附近的一個(gè)巨大渦，舵的上下表面幾乎都沒有渦形成。由于已經(jīng)忽略了阻尼的作用，所以舵的俯仰運(yùn)動(dòng)繼續(xù)加劇，而非經(jīng)典的等幅運(yùn)動(dòng)。這不僅會(huì)大大降低航行體的操控性能，甚至有可能造成舵系統(tǒng)的結(jié)構(gòu)損壞。

圖12 V∞=2.37 m/s時(shí)，舵附近流場(chǎng)以及渦量云圖Fig.12 Vorticity contours of the hydrofoil at V∞=2.37 m/s

3 結(jié) 論

本文通過對(duì)經(jīng)典兩自由度系統(tǒng)任意運(yùn)動(dòng)的無量綱動(dòng)力方程進(jìn)行參數(shù)修正，在計(jì)算中引入了展弦比對(duì)流體動(dòng)力的影響。兩自由度舵系統(tǒng)的流致振動(dòng)實(shí)際上是動(dòng)能和勢(shì)能相互轉(zhuǎn)化的過程，隨著V∞的增加相軌線的形狀越來越趨于扁平化，這意味著系統(tǒng)的動(dòng)能可以轉(zhuǎn)化為更大的勢(shì)能。這是因?yàn)榈土魉贂r(shí)，系統(tǒng)的動(dòng)能會(huì)產(chǎn)生耗散，致使系統(tǒng)進(jìn)行衰減振動(dòng)；當(dāng)流速達(dá)到VF時(shí)，系統(tǒng)的動(dòng)能可以全部轉(zhuǎn)化為勢(shì)能，勢(shì)能亦可全部轉(zhuǎn)化為動(dòng)能，從而出現(xiàn)極限環(huán)形式的相軌線；而在高流速下，系統(tǒng)的動(dòng)能不但不會(huì)發(fā)生耗散，還會(huì)在振動(dòng)過程中不斷從流體中吸收能量并將其轉(zhuǎn)化為勢(shì)能，使得系統(tǒng)振動(dòng)相軌線不斷向外擴(kuò)散。

低質(zhì)量比舵系統(tǒng)發(fā)生顫振時(shí)，系統(tǒng)的振動(dòng)形態(tài)從穩(wěn)態(tài)轉(zhuǎn)化為非穩(wěn)態(tài)。在對(duì)應(yīng)的頻譜圖上，系統(tǒng)的沉浮頻率或俯仰頻率附近會(huì)出現(xiàn)一個(gè)新的頻率，表明在這個(gè)過渡期間系統(tǒng)的運(yùn)動(dòng)受到了一定程度的擾動(dòng)，但是并沒有像高質(zhì)量比系統(tǒng)那樣形成一個(gè)突出的頻率，這是低質(zhì)量比系統(tǒng)和高質(zhì)量比系統(tǒng)顫振現(xiàn)象的主要區(qū)別。

通過數(shù)值模擬，發(fā)現(xiàn)在舵系統(tǒng)運(yùn)動(dòng)最初，舵的后緣附近會(huì)出現(xiàn)一系列細(xì)小的尾渦。并且，在舵的前緣會(huì)出現(xiàn)了明顯的流體分離，這種分離不斷向后緣擴(kuò)散，同時(shí)會(huì)在舵的尾部出現(xiàn)一系列脫落渦。當(dāng)系統(tǒng)進(jìn)行衰減運(yùn)動(dòng)時(shí)，則尾渦會(huì)隨著系統(tǒng)的靜止而消失，舵上表面的渦層主要出現(xiàn)在1/4弦長(zhǎng)處和后緣之間。如果系統(tǒng)進(jìn)行發(fā)散運(yùn)動(dòng)，則舵會(huì)發(fā)生劇烈運(yùn)動(dòng)，其前緣附近產(chǎn)生巨大的渦，渦層遍布舵的整個(gè)上表面。并且在1/4弦長(zhǎng)處會(huì)出現(xiàn)明顯層流分離現(xiàn)象，這與衰減運(yùn)動(dòng)不同。