兩自由度運(yùn)動(dòng)圓柱繞流的離散渦方法模擬

董 婧，宗 智，李章銳，孫 雷，陳 偉

（1大連理工大學(xué) 工業(yè)裝備結(jié)構(gòu)分析國家重點(diǎn)實(shí)驗(yàn)室 船舶與海洋工程學(xué)院，遼寧 大連116024；

2中國艦船研究設(shè)計(jì)中心，武漢430050）

1 引 言

將鈍體放置在具有一定速度的來流當(dāng)中，流體會(huì)繞過結(jié)構(gòu)物，在其后面形成尾流區(qū),并在一定條件下從結(jié)構(gòu)兩側(cè)交替地發(fā)生周期性泄渦現(xiàn)象。而旋渦的產(chǎn)生和脫落的作用，使物體受到橫向和流向的脈動(dòng)壓力。對于彈性圓柱體或彈性支承的剛性圓柱體，所產(chǎn)生的脈動(dòng)力會(huì)引起結(jié)構(gòu)物的振動(dòng)，反過來結(jié)構(gòu)物的振動(dòng)也會(huì)影響到流場，這種流固耦合現(xiàn)象稱為渦激振動(dòng)（Vortex Induced Vibration）[1]。當(dāng)渦泄頻率接近于結(jié)構(gòu)的某一固有頻率時(shí)，旋渦的泄放過程被結(jié)構(gòu)的振動(dòng)所控制,而使旋渦泄放和結(jié)構(gòu)振動(dòng)具有相同的頻率，并在較大的流速范圍內(nèi)保持不變，這種現(xiàn)象稱為“鎖定”[2-3]。在發(fā)生鎖定時(shí)物體的振幅顯著增加，相位角發(fā)生突變。鎖定現(xiàn)象會(huì)在旋渦脫落頻率的很大范圍內(nèi)發(fā)生，鎖定現(xiàn)象的機(jī)理尚不明確，大量的實(shí)驗(yàn)數(shù)據(jù)表明，通常當(dāng)折合速度（Ur）在4.5到8之間時(shí)發(fā)生鎖定現(xiàn)象。

海洋工程中的水下結(jié)構(gòu)物經(jīng)常采用圓柱體，其在波流作用下經(jīng)常會(huì)發(fā)生渦激振動(dòng)。當(dāng)泄渦頻率接近于結(jié)構(gòu)的某一固有頻率時(shí)，將會(huì)出現(xiàn)垂直于水流和圓柱體軸線方向的強(qiáng)烈振動(dòng)，并且導(dǎo)致沿水流方向的振動(dòng)和曳力增大。結(jié)構(gòu)在兩個(gè)方向上的強(qiáng)烈振動(dòng)將引起震蕩應(yīng)力[4]。如果這種振動(dòng)處于定常狀態(tài)而連續(xù)發(fā)生，將會(huì)導(dǎo)致結(jié)構(gòu)的疲勞破壞，降低其使用壽命。渦激振動(dòng)是海洋立管遭受破壞的主要原因，另外，在水流沖刷及海底表面不平等環(huán)境因素作用下，鋪設(shè)在海底的油氣管道難免形成懸空段，直接暴露在水流作用下也會(huì)發(fā)生渦激振動(dòng)。

關(guān)于結(jié)構(gòu)渦激振動(dòng)的研究，國外起于20世紀(jì)60年代，而我國則從80年代開始。研究范圍涉及風(fēng)工程和海洋工程；研究方法在早期側(cè)重于試驗(yàn)研究，經(jīng)過試驗(yàn)研究和理論分析并重的階段，目前有重視理論分析的傾向[2-4]。眾所周知，在通常的渦激振動(dòng)中橫向振動(dòng)占主導(dǎo)地位，振幅大大超過流向振動(dòng)，所以在研究過程中，很多學(xué)者采用簡化的單自由度運(yùn)動(dòng)模型，只考慮橫向振動(dòng)而忽略了流向振動(dòng)因素。那么這種做法是否科學(xué)合理？流向振動(dòng)到底對渦激振動(dòng)有何影響？本文的目的之一就要驗(yàn)證它。此外，本文還要對渦激振動(dòng)的特征機(jī)理進(jìn)行分析探討。

本文采用離散渦方法（Discrete Vortex Method，DVM）來研究彈性支承的二維圓柱繞流的渦激振動(dòng)（VIV）問題。對于圓柱的運(yùn)動(dòng)情況的研究分為單向自由度的橫向運(yùn)動(dòng)和兩向自由度的流向和橫向共同運(yùn)動(dòng)兩種情況。本文分別進(jìn)行了單向自由度橫向運(yùn)動(dòng)和兩向自由度橫向和流向共同運(yùn)動(dòng)這兩種情況的模擬，并討論了兩者之間的聯(lián)系，在一系列共16組（單向、兩向運(yùn)動(dòng)情況各8組）數(shù)值計(jì)算中采用同Stappenbelt[6]的試驗(yàn)一致的各項(xiàng)參數(shù)，其中包括相同尺度的圓柱模型，均為單位長度，直徑0.055 4m；從2.36到12.96變化的8個(gè)不同的質(zhì)量比，以及隨質(zhì)量比的改變而改變的水中固有頻率；從3到16變化的折合速度，以及隨折合速度和固有頻率的改變而改變的速度和雷諾數(shù)；還有不變的阻尼比，按試驗(yàn)值均為0.006。

本文首先對受力系數(shù)和響應(yīng)位移進(jìn)行了時(shí)域分析。從脈動(dòng)升力和圓柱橫向響應(yīng)的振幅、周期和相位關(guān)系來判斷鎖定區(qū)域和非鎖定區(qū)域，討論在不同折合速度下，振動(dòng)幅值的初始分支、上端分支和下端分支的變化規(guī)律，以及質(zhì)量比的增加對振動(dòng)幅值、鎖定位置和鎖定區(qū)域的影響。

接下來進(jìn)行的是頻域響應(yīng)的分析，就泄渦頻率和振動(dòng)頻率關(guān)系與Stappenbelt實(shí)驗(yàn)結(jié)果進(jìn)行了對比，吻合完好，斯托哈爾數(shù)都處于穩(wěn)定范圍內(nèi)，并且進(jìn)一步討論了質(zhì)量比的增加對振動(dòng)幅值、鎖定位置和鎖定區(qū)域的影響。

最后進(jìn)行了不同折合速度、不同自由度下尾渦模式的分析，并結(jié)合一些文獻(xiàn)中已發(fā)表的結(jié)果進(jìn)行了對比和討論。在分析的過程中認(rèn)為數(shù)值模擬的絕大部分尾渦模式與實(shí)驗(yàn)和計(jì)算相一致，并且發(fā)現(xiàn)了新的尾渦模式。

2 離散渦方法的原理[1，5]

假設(shè)流體為不可壓縮流體，那么它的控制方程為N-S方程：

因渦量ω=▽×V，對上式取旋度，得：

如只研究二維平面流，則ω·▽V=0，由此得：

離散渦方法是將連續(xù)的渦量場離散，在拉格朗日框架下進(jìn)行求解，著眼于追蹤每個(gè)渦元粒子運(yùn)動(dòng)的時(shí)間歷程。對于不可壓縮的粘性流動(dòng)，上世紀(jì)90年代Chorin[7]提出了算子分裂法，渦量方程（2）可分裂成對流和擴(kuò)散兩部分分別進(jìn)行求解。

對流部分：

粘性擴(kuò)散部分：

其中，υ為流體的運(yùn)動(dòng)粘性系數(shù)，ω為渦量，ψ為流函數(shù)。

方程（4）代表二維不可壓縮無粘性流動(dòng)，如果當(dāng)無窮遠(yuǎn)處的自由流動(dòng)條件已知，流場中不存在外部邊界時(shí)，它的解就是著名的Biot-Savart定律：

式中 r、r′，表示矢徑，V∞代表無窮遠(yuǎn)處的速度。

而擴(kuò)散部分與粒子的布朗運(yùn)動(dòng)的解具有相同的形式，在無粘解的基礎(chǔ)上，可以采用隨機(jī)走位的方法來模擬，下面具體介紹隨機(jī)走位的方法。

方程（5）為一維熱擴(kuò)散型方程，它的基本解是格林函數(shù)：

這與一個(gè)均方差為零、標(biāo)準(zhǔn)差為σx的隨機(jī)變量的概率密度函數(shù)在形式上是一致的，即：

其中N是流場中渦元總數(shù)，下標(biāo)i為渦元編號(hào)，Pi和Qi分別是在（0,1）和（0，2π）區(qū)間內(nèi)均勻分布的兩個(gè)相互獨(dú)立的隨機(jī)數(shù)。

假設(shè)t時(shí)刻渦元的位置xi（t），yi（t）已給定，則t+Δt時(shí)刻渦元的隨機(jī)走位為

在此應(yīng)當(dāng)指出的是，每個(gè)渦元的環(huán)量Γi是不隨時(shí)間變化的。

3 離散渦方法的數(shù)值實(shí)現(xiàn)

離散渦方法的實(shí)質(zhì)是用離散的點(diǎn)渦來模擬物體的邊界層和剪切層，認(rèn)為粘性影響集中在物體附近很薄的邊界層以及分離后的剪切層內(nèi)，其渦量來自于不斷分離出來的點(diǎn)渦供給，在整個(gè)流場中渦量守恒。離散渦方法主要分為物面渦的生成、渦元的對流和擴(kuò)散、渦元的消除和融合，物體受力計(jì)算等幾個(gè)方面。

首先，將連續(xù)的渦量離散成為若干個(gè)渦元（如圖1所示），根據(jù)流場條件計(jì)算出每個(gè)渦元的位置和強(qiáng)度。

其次，要求出流場中的速度分布，每個(gè)渦元的速度都由兩部分疊加而成，一部分為定常流場的繞流速度；另一部分為渦元之間的誘導(dǎo)速度，可以通過泊松方程▽2ψ=-ω和邊界條件聯(lián)立求解。

在求解邊界條件時(shí)，本文采用流函數(shù)列式的方法即渦元法（Vortex Element Method）[1，5]。 渦元法利用物面是一條流線

圖1 離散渦元示意圖Fig.1 Position of vortex generation

來滿足物面條件，通過離散的流函數(shù)列式，可以確定各新生渦元的強(qiáng)度。Spalart等人在1983年引入了剛性渦，其半徑為σ，周向的誘導(dǎo)速度為：

當(dāng)渦元的誘導(dǎo)速度為（5）式時(shí)，其相應(yīng)的流函數(shù)為

則在控制點(diǎn)k處的流函數(shù)為：

其中N和M分別為新生渦元總數(shù)和尾渦渦元總數(shù)。

當(dāng)圓柱靜止不動(dòng)時(shí)物面為一條流線，則

將（13）式代入（14）式可以得到一組以環(huán)量Γ為變量的線性方程組：

其中系數(shù)矩陣AI僅與控制點(diǎn)和新生渦元的位置有關(guān)，列向量b由來流和尾流中的渦元分布來決定。

關(guān)于物體的受力，Quartepelle和Napolitanp在1983年推導(dǎo)了用整個(gè)流場的信息來表示的受力系數(shù)公式[9]：

對于運(yùn)動(dòng)的圓柱，（14）式應(yīng)改寫為：

其中VB圓柱運(yùn)動(dòng)的速度矢量，n為控制點(diǎn)k處的物面法向。相應(yīng)的（15）式中的列向量b中也會(huì)增加與速度有關(guān)的因素，而系數(shù)矩陣AI保持不變。此時(shí)（16）式必須加入圓柱的慣性力項(xiàng)，于是改寫為：

無論圓柱運(yùn)動(dòng)與否，現(xiàn)有渦元的速度都由環(huán)量Γ表示如下：從上式中可以看出，在求取每個(gè)渦元的速度時(shí)，都需要計(jì)算其他所有渦元對它的誘導(dǎo)速度，如果把計(jì)算單個(gè)渦元速度的運(yùn)算次數(shù)計(jì)為NV，那么在一個(gè)時(shí)間步內(nèi)計(jì)算所有渦元的速度，則需要N2V次。而每一時(shí)間步內(nèi)都有固定數(shù)目的新生渦元進(jìn)入尾流中，隨著時(shí)間的增長NV直線增加，其結(jié)果必然導(dǎo)致單個(gè)時(shí)間步長的計(jì)算時(shí)間越來越長，計(jì)算緩慢。因此，必須有效控制尾流中渦元的數(shù)目。同時(shí)，在渦的對流和擴(kuò)散中，渦元會(huì)越過邊界移動(dòng)到圓柱內(nèi)部引起混亂，也必須予以消除。Spalart提出的渦融合機(jī)制認(rèn)為，如果兩渦元滿足下面（20）式的條件就進(jìn)行合并。

其中z=x+yi為渦元在復(fù)平面內(nèi)的位置，D0和V0都是控制參數(shù)，其中D0∝U∞，控制圓柱附近的渦元數(shù)目；V0∝10-6U∞，控制流場中渦元總數(shù)。合并后新渦元的位置為（21）式，強(qiáng)度為（22）式。

運(yùn)動(dòng)圓柱的渦激振動(dòng)問題可以從兩方面入手，即強(qiáng)迫振動(dòng)和自激振動(dòng)。本文進(jìn)行了自激振動(dòng)的模擬，分別采用單向自由度和兩向自由度的彈簧-阻尼質(zhì)量系統(tǒng)來模擬二維圓柱，如圖2所示。

圓柱在橫向和流向的運(yùn)動(dòng)方程可分別表示為：

圖2 彈性支承的圓柱模型Fig.2 Vibration model of elastic mounted cylinder

其中m為圓柱的質(zhì)量，c為系統(tǒng)阻尼系數(shù)，k為彈簧剛度系數(shù)。Fx（t），F(xiàn)y（t）分別為流向和橫向的流體力，它們可以分別用阻力系數(shù)和升力系數(shù)表示為：

將（23）、（24）式轉(zhuǎn)化成無量綱形式為：

采用4階精度的Runge-Kutta法進(jìn)行迭代求解。

4 數(shù)值結(jié)果和實(shí)驗(yàn)對比

在兩向自由度橫向振動(dòng)的數(shù)值模擬中，為方便與其它試驗(yàn)和數(shù)值結(jié)果對比，本文選用Stappenbelt[3]的實(shí)驗(yàn)參數(shù)，如表1所示

除表1所列的參數(shù)之外還有其他一些在這8組實(shí)驗(yàn)當(dāng)中不變的參數(shù)，包括圓柱直徑D為0.055 4m，折合速度Ur從3變化到16，以0.5遞增，時(shí)間步長Δt=0.1D/U，每時(shí)間步有新生渦120個(gè)，總共1 200步。雷諾數(shù)Re隨折合速度Ur變化，其值在2.84×105到1.52×106之間，每組實(shí)驗(yàn)都分為單自由度橫向振動(dòng)和兩自由度流向和橫向耦合振動(dòng)兩種情況計(jì)算，為方便表述，分別記為1dof和2dof。

在圖 3（a）到（f）給出了質(zhì)量比為 2.36時(shí)的 6個(gè)典型時(shí)刻的升力系數(shù)CL、阻力系數(shù)CD、無因次化的橫向位移Y/D、流向位移X/D的時(shí)程曲線。通觀這6幅圖可以得到以下結(jié)論：

表1 振動(dòng)模型計(jì)算參數(shù)Tab.1 Parameter of vibration model

從受力系數(shù)來看，阻力系數(shù)CD始終處于一個(gè)較高的水平，其波動(dòng)的幅值較小，這其中穩(wěn)定的部分是由均勻來流的持續(xù)沖擊所致，而波動(dòng)的部分是由泄渦產(chǎn)生的脈動(dòng)壓力所致，所以從上面的幾幅圖的比較可以看出，折合速度越大，阻力系數(shù)的平均值就越高，這是符合常理的。升力系數(shù)CL的平衡位置在0附近，而波動(dòng)的幅值相對較大，這是由于升力系數(shù)只與泄渦產(chǎn)生的脈動(dòng)壓力有關(guān)，橫向的脈動(dòng)壓力幅值大于流向的脈動(dòng)壓力幅值，而均勻來流的持續(xù)沖擊力對升力系數(shù)沒有直接的影響。升力系數(shù)CL的脈動(dòng)周期始終約為阻力系數(shù)CD脈動(dòng)周期的2倍，這也是由泄渦產(chǎn)生的脈動(dòng)壓力所致。從位移曲線來看，流向位移X/D受阻力系數(shù)CD影響，在一個(gè)較高的水平上波動(dòng)，在第一個(gè)周期里有較高的幅值，這是由于在靜止?fàn)顟B(tài)突然加力所致，在后面的幾個(gè)周期里就會(huì)平穩(wěn)下來，不過這也在另一個(gè)側(cè)面說明了本文計(jì)算程序的穩(wěn)定性。橫向位移Y/D受升力系數(shù)CL影響，則在0軸上下波動(dòng)。

圖3（a）中曲線對應(yīng)折合速度Ur=4.5，處于初始分支。此時(shí)升力系數(shù)CL和阻力系數(shù)CD出現(xiàn)階段性的脈動(dòng)，周期大概在25左右，這與彈簧—阻尼系統(tǒng)在兩個(gè)方向的固有頻率接近有關(guān)。橫向位移Y/D與升力系數(shù)CL的周期和相位都是一致的，這一現(xiàn)象在初始分支中普遍存在，并與眾多試驗(yàn)和數(shù)值模擬結(jié)論相一致[10-12]。圖3（b）中曲線對應(yīng)折合速度Ur=6，處于初始分支和上端分支的臨界狀態(tài)，橫向位移Y/D的振動(dòng)幅值較大且振動(dòng)平穩(wěn)，但相位和周期仍然符合初始分支的特點(diǎn)。圖3（c）中曲線對應(yīng)折合速度Ur=7.5，處于上端分支，橫向位移Y/D的振動(dòng)幅值較大且振動(dòng)平穩(wěn)，橫向位移Y/D與升力系數(shù)CL的周期一致而相位角相差180°，這是上端分支鎖定區(qū)域的典型特點(diǎn)。圖3（d）中曲線對應(yīng)折合速度Ur=9，處于上端分支和下端分支的臨界狀態(tài)，此時(shí)橫向位移Y/D與升力系數(shù)CL仍然是周期一致而相位角相差180°，而橫向位移Y/D的振動(dòng)幅值開始下降且振動(dòng)不穩(wěn)定。圖3（e）、（f）處于下端分支的中段和末段，在這一階段，橫向振動(dòng)逐漸減退，流向振動(dòng)和橫向振動(dòng)的幅值相當(dāng)。

圖3 升力系數(shù)CL、阻力系數(shù)CD和流向位移X/D、橫向位移Y/D的時(shí)程曲線，其中質(zhì)量比m*=2.36，水中固有頻率Fn=1.711Fig.3 History of lift coefficient CL,drag coefficient CD,stream-wise displacement X/D,and transverse displacement Y/D,here mass-ratio m*=2.36,natural frequency in water Fn=1.711

圖4到圖7為不同質(zhì)量比時(shí)單自由度1dof和兩自由度2dof的橫向運(yùn)動(dòng)無因次振幅A/D同實(shí)驗(yàn)值的對比，橫坐標(biāo)為折合速度Ur。為與實(shí)驗(yàn)的描述一致，這里用T表示橫向運(yùn)動(dòng)，T&S表示橫向和流向耦合運(yùn)動(dòng)。從圖中可以看出，本文計(jì)算結(jié)果與Stappenbelt實(shí)驗(yàn)值的總體趨勢向一致，會(huì)出現(xiàn)明顯的初始分支、上端分支和下端分支。在兩自由度的振動(dòng)中質(zhì)量比較低時(shí)還會(huì)出現(xiàn)超上端分支，而后振幅驟降。從相同質(zhì)量比時(shí)1dof和2dof振動(dòng)情況的比較來看，相同折合速度下2dof振動(dòng)時(shí)的幅值更大，這一特征隨著質(zhì)量比的減小而增加，這說明了流向振動(dòng)的參與對橫向振動(dòng)會(huì)有促進(jìn)作用。而從不同質(zhì)量比之間的對比可以看出，隨著質(zhì)量比的增加，無論是1dof還是2dof振動(dòng)，其幅值都會(huì)隨之下降。本文計(jì)算結(jié)果的不足之處在于模擬到的下端分支處的振幅偏高。

圖4 m*=2.36時(shí)單雙自由度橫向無因次振幅A/D同實(shí)驗(yàn)對比Fig.4 Non-dimensional transverse amplitudes A/D of 1dof and 2dof compared with experimental results for m*=2.36

圖5 m*=3.68時(shí)單雙自由度橫向無因次振幅A/D同實(shí)驗(yàn)對比Fig.5 Non-dimensional transverse amplitudes A/D of 1dof and 2dof compared with experimental results for m*=3.68

圖6 m*=5.19時(shí)單雙自由度橫向無因次振幅A/D同實(shí)驗(yàn)對比Fig.6 Non-dimensional transverse amplitudes A/D of 1dof and 2dof compared with experimental results for m*=5.19

圖7 m*=6.54時(shí)單雙自由度橫向無因次振幅A/D同實(shí)驗(yàn)對比Fig.7 Non-dimensional transverse amplitudes A/D of 1dof and 2dof compared with experimental results for m*=6.54

圖8共有10組質(zhì)量比的單自由度運(yùn)動(dòng)時(shí)無因次振動(dòng)頻率和泄渦頻率圖，橫坐標(biāo)為折合速度。其中前兩幅圖為Stappenbelt實(shí)驗(yàn)中已發(fā)表的數(shù)據(jù)。后8幅圖為本文計(jì)算結(jié)果，與實(shí)驗(yàn)結(jié)果的趨勢一致，在鎖定區(qū)域以內(nèi)，泄渦頻率fs和圓柱振動(dòng)頻率fv一起鎖定在固有頻率fn附近；在鎖定區(qū)域以外，泄渦頻率fs沿折合速度會(huì)形成一條直線，在這條直線上有一個(gè)相對穩(wěn)定的斯托哈爾數(shù)St，這將在表2中說明。

表2為單自由度振動(dòng)模型的鎖定參數(shù) （包括初始鎖定的折合速度Ur和鎖定區(qū)域折合速度的范圍ΔUr）及斯托哈爾數(shù)St的計(jì)算結(jié)果同Stappenbelt實(shí)驗(yàn)值的對比。從折合速度Ur的比較中可以看出，兩組數(shù)據(jù)的大小相接近，且都隨著質(zhì)量比的增加而增加。從鎖定區(qū)域的范圍ΔUr可以看出，本文計(jì)算與實(shí)驗(yàn)結(jié)果符合得較好，最多相差0.6，且都隨著質(zhì)量比的增加而減小。從斯托哈爾數(shù)St的比較可以看出，本文計(jì)算與實(shí)驗(yàn)結(jié)果符合得較好，最多相差0.08，且都隨著質(zhì)量比的增加而增大。雙自由度振動(dòng)模型的計(jì)算結(jié)果與之類似。

圖8 無因次振動(dòng)頻率和泄渦頻率關(guān)系圖Fig.8 Non-dimensional vibration and vortex shedding frequencies

表2 鎖定及斯托哈爾數(shù)計(jì)算結(jié)果同實(shí)驗(yàn)對比Tab.2 Lock-in values and Strouhal number compared with experimental results

橫向振動(dòng)響應(yīng)分支在振幅上有較大變化，相應(yīng)的尾渦模式也有明顯的差異。在本文的計(jì)算結(jié)果當(dāng)中，Re在105以上，無論是1dof還是2dof的振動(dòng)，在初始分支中都表現(xiàn)出2S模式（two single vortex shedding per cycle）如圖9所示，而通常情況下在上端分支和下端分支所對應(yīng)的尾渦則為2P模式（Two pairs vortex shedding per cycle）。上端分支所對應(yīng)的尾渦中2P模式為一強(qiáng)一弱的兩個(gè)旋渦如圖10所示，其中小的旋渦會(huì)在尾流中不斷的對流和擴(kuò)散中消失，在流場中較遠(yuǎn)處只剩下較大的旋渦；而下端分支所對應(yīng)的尾渦中2P模式為強(qiáng)度相當(dāng)?shù)膬蓚€(gè)旋渦，如圖11所示，這兩個(gè)漩渦在流場較大范圍內(nèi)運(yùn)動(dòng)都不會(huì)消失。這與黃志勇[11]文獻(xiàn)中描述的結(jié)果類似。而本文的計(jì)算中還發(fā)現(xiàn)了，在個(gè)別上端分支中在一個(gè)周期內(nèi)先是產(chǎn)生兩條細(xì)長的渦帶而后分化成4個(gè)小的旋渦，暫時(shí)把這種尾渦模式也歸為2P模式的一種，如圖12所示，本文認(rèn)為這是從2P模式到2T模式的一種過渡模式，其生成機(jī)理有待于進(jìn)一步研究。

圖9 初始分支的2S瀉渦模式Fig.9 2S model of vortex shedding in the initial branch

圖10 上端分支中兩旋渦不等的2P瀉渦模式Fig.10 2P model of vortex shedding in the upper branch(each pair of vortex blobs with unequal scales)

圖11 下端分支中兩旋渦相等的2P瀉渦模式Fig.11 2P model of vortex shedding in the lower branch(each pair of vortex blobs with equal scales)

圖12 特殊的2P瀉渦模式Fig.12 The special 2P model of vortex shedding

圖13 超上端分支中的2T瀉渦模式Fig.13 2T model of vortex shedding in the super upper branch

2T模式（Two triplets of vortices per cycle）在一個(gè)泄渦周期內(nèi)產(chǎn)生兩組渦團(tuán)，每個(gè)渦團(tuán)產(chǎn)生三個(gè)旋渦，其中兩個(gè)旋渦的旋向與另一個(gè)相反。出現(xiàn)在兩自由度振動(dòng)的上端分支的頂點(diǎn)處，也就是所說的超上端分支，如圖13所示。

在一些過渡階段還會(huì)出現(xiàn)以上多種尾渦模態(tài)交替出現(xiàn)的情況，是由于運(yùn)動(dòng)圓柱的位移和速度的變化所致，在這里不做詳細(xì)介紹。

5 結(jié) 論

本文利用離散渦方法對二維彈性支承圓柱繞流渦激振動(dòng)問題進(jìn)行了較為詳細(xì)的計(jì)算，可以得到以下結(jié)論：

首先對受力系數(shù)和響應(yīng)位移進(jìn)行了時(shí)域分析。從脈動(dòng)升力和圓柱橫向響應(yīng)的周期和相位關(guān)系可以看出有明顯的鎖定區(qū)域和非鎖定區(qū)域，觀察到在不同折合速度下，振動(dòng)幅值有著明顯的初始分支、上端分支和下端分支的變化，并且隨著質(zhì)量比的增加，振動(dòng)幅值會(huì)有所下降，這與眾多觀察和實(shí)驗(yàn)的本質(zhì)特征一致。

在單自由度和兩自由度振動(dòng)幅值的對比中發(fā)現(xiàn)，兩自由度振動(dòng)在上端分支處的幅值明顯高于單自由度的情況，說明流向振動(dòng)對橫向振動(dòng)有促進(jìn)作用，用單自由度建立的渦激振動(dòng)模型預(yù)報(bào)結(jié)果將會(huì)偏于保守，在此建議在建立渦激振動(dòng)模型時(shí)，流向振動(dòng)因素應(yīng)當(dāng)予以重視。

從進(jìn)一步的頻域分析結(jié)果來看，泄渦頻率和振動(dòng)頻率關(guān)系與Stappenbelt實(shí)驗(yàn)結(jié)果吻合完好，斯托哈爾數(shù)都處于穩(wěn)定范圍內(nèi)，并且驗(yàn)證了隨著質(zhì)量比的增加，初始鎖定時(shí)的折合速度會(huì)提高，而鎖定區(qū)域的范圍會(huì)減小。

在進(jìn)行尾渦模式分析時(shí)與黃志勇在文獻(xiàn)中已發(fā)表的結(jié)論基本一致，并且在上端分支中發(fā)現(xiàn)了新的尾渦模式—有4個(gè)小的旋渦的2P模式。本文認(rèn)為這是從2P模式到2T模式的一種過渡模式。

[1]Lewis R I.Vortex element method for fluid dynamic analysis of engineering systems[M].Cambridge:Cambridge University Press,1991.

[2]Gabbai R D,Benaroya H.An overview of modeling and experiments of vortex-induced vibration of circular cylinders[J].Journal of Sound and Vibration,2005,282(3):575-616.

[3]Liao Jungchi.Vortex-induced vibration of slender structures in unsteady flow[M].Cambridge,USA:Massachusetts Institute of Technology,2002.

[4]Pan Zhiyuan,Cui Weicheng,Liu Yingzhong.A predicting model for self-excited VIV of a circular cylinder at low massdamping[J].Journal of Ship Mechanics,2005,10(5):115-124.(in Chinese)

[5]Saltara F.Meneghini J R,Fregonesi R A.Numerical simulation of flow around elastically mounted cylinder[J].International Journal of Offshore and Polar Engineering,2003,13(2):99-104.

[6]Stappenbelt B,Lalji F.Low mass ratio vortex-induced motion[C].16th Australasian Fluid Mechanics Conference,2007,12:1491-1497.

[7]Chorin A J.Numerical study of slightly viscous flow[J].Journal of Fluid Mechanics,1973,57:785-796.

[8]Lewis R I.Vortex element method for fluid dynamic analysis of engineering systems[M].Cambridge:Cambridge University Press,1991.

[9]Chorin A J.Numerical study of slightly viscous flow[J].Journal of Fluid Mechanics,1973,57:785-796.

[10]Huang Zhiyong,Pan Zhiyuan,Cui Weicheng.Numerical simulation of VIV of a circular cylinder with two degrees of freedom and low mass-ratio[J].Journal of Ship Mechanics,2007,11(1):1-9.(in Chinese)

[11]Khalak A,Williamson C H K.Dynamics of a hydroelastic cylinder with very low mass and damping[J].Journal of Fluids and Structures,1999,13:813-851.

[12]Zhou C Y,So R M C,Lam K.Vortex-induced vibrations of an elastic circular cylinder[J].Journal of Fluids and Structures,1999,13:165-189.