及春寧, 譚沛森， 仵柏衡

(1.天津大學(xué)水利工程仿真與安全國(guó)家重點(diǎn)實(shí)驗(yàn)室，天津 300072； 2.中國(guó)海洋大學(xué)海洋與大氣學(xué)院，山東 青島 266100； 3.挪威科技大學(xué)海洋技術(shù)系，挪威 特隆赫姆 7491)

振蕩流中的靜止圓柱以及靜止流中的振蕩圓柱在海洋工程和水動(dòng)力學(xué)領(lǐng)域有著極大研究?jī)r(jià)值。該物理模型適用于很多海洋工程設(shè)施，如研究液化天然氣運(yùn)輸船在受到波浪力擺動(dòng)時(shí)，液艙中泵塔的受力；由多個(gè)圓柱組成的波能吸收器在波浪中的振動(dòng)；水下防噴器受到擾動(dòng)在流速較小的海水中的振蕩[1]等。對(duì)于上述工程應(yīng)用，渦激振動(dòng)產(chǎn)生的力學(xué)問(wèn)題會(huì)對(duì)建筑設(shè)施產(chǎn)生極大的影響，忽視這種影響可能釀成災(zāi)害性后果，因此研究圓柱體在振蕩過(guò)程中的受力載荷至關(guān)重要。實(shí)際工程中，圓柱的數(shù)量和擺放方式各異，而并列排布的雙圓柱系統(tǒng)是其中較為典型的一種。對(duì)雙圓柱體在振蕩過(guò)程中的受力載荷研究對(duì)于大量的實(shí)際工程應(yīng)用問(wèn)題具有突出的指導(dǎo)意義。

該問(wèn)題的流體動(dòng)力學(xué)模型可以簡(jiǎn)化為振蕩流體中的結(jié)構(gòu)物相互作用。該結(jié)構(gòu)上沿振動(dòng)方向的力可以由莫里森方程表示:

而垂直于振動(dòng)方向的力可以表示為：

其中：Cd為阻力系數(shù)；Cm為附加質(zhì)量系數(shù)；Cl為升力系數(shù)。Keulegan[11]指出，這三個(gè)關(guān)鍵的水動(dòng)力學(xué)系數(shù)依賴于無(wú)量綱數(shù)KC。KC數(shù)被定義為KC=UmT/D，Um為圓柱最大的振蕩速度，T是圓柱的振蕩周期，D為圓柱直徑。在簡(jiǎn)單正弦振動(dòng)中，KC數(shù)可以簡(jiǎn)化為2πA/D,A為圓柱振蕩的振幅。Sarpkaya[12]指出，這三個(gè)水動(dòng)力學(xué)系數(shù)不僅依賴于KC數(shù)，也同時(shí)依賴于另一個(gè)無(wú)量綱數(shù)β，β被定義為β=D2/Tν，ν為流體的運(yùn)動(dòng)學(xué)粘度系數(shù)。因?yàn)槔字Z數(shù)Re定義為Re=UmD/v因此雷諾數(shù)Re也可表示為Re=β×KC。

單個(gè)圓柱在均勻流中的水動(dòng)力特性已經(jīng)得到了廣泛的研究。Bearman等[3]實(shí)驗(yàn)測(cè)量了KC數(shù)范圍為4～55和β范圍為100～166時(shí)的圓柱體上的流體力合力，小KC數(shù)范圍內(nèi)，阻力系數(shù)Cd與理論結(jié)果吻合較好。Sarpkaya[4]實(shí)驗(yàn)在穩(wěn)定流場(chǎng)所對(duì)應(yīng)的KC數(shù)范圍下，合理的預(yù)測(cè)阻力系數(shù)和附加質(zhì)量系數(shù)隨著KC數(shù)上升的變化趨勢(shì)。Tatsuno和Bearman[5]在不同的KC數(shù)值下進(jìn)行了實(shí)驗(yàn)，并確定了8種不同的流動(dòng)狀態(tài)。

與單個(gè)圓柱相比，多體圓柱在靜止流場(chǎng)中振蕩的水動(dòng)力特性會(huì)復(fù)雜很多。這是圓柱的尾流之間的相互作用所導(dǎo)致的。對(duì)于串聯(lián)圓柱系統(tǒng)，下游圓柱受到上游圓柱尾流作用，在間隙比較小時(shí)會(huì)產(chǎn)生馳振現(xiàn)象(Wake-induced galloping)[6]。Chen等[7]對(duì)雷諾數(shù)為100，間隙比為2.0～5.0的雙圓柱系統(tǒng)進(jìn)行了數(shù)值模擬，并分析了不對(duì)稱振動(dòng)的現(xiàn)象。為了進(jìn)一步研究雙圓柱系統(tǒng)水動(dòng)力特征和其受力載荷，本文首先進(jìn)行了了靜水中振蕩的并列雙圓柱動(dòng)力特性的模型實(shí)驗(yàn)。此外，為了更好的理解雙圓柱的流場(chǎng)特征及其機(jī)理，本文同時(shí)對(duì)并列雙圓柱系統(tǒng)在靜止流場(chǎng)中的運(yùn)動(dòng)進(jìn)行了數(shù)值模擬，從能量的角度對(duì)阻力系數(shù)增大進(jìn)行了機(jī)理分析。

1 實(shí)驗(yàn)設(shè)置

1.1 實(shí)驗(yàn)裝置與器材

本研究實(shí)驗(yàn)在麻省理工學(xué)院拖曳水池實(shí)驗(yàn)室的拖曳水池中進(jìn)行，研究靜水中振蕩圓柱的水動(dòng)力特性，實(shí)驗(yàn)整體布置如圖1(a)所示,水池長(zhǎng)1.5 m，寬和高均為1 m。實(shí)驗(yàn)所使用圓柱的表面經(jīng)過(guò)電鍍，以保持其表面光滑，從而減小因表面粗糙引起的誤差，圓柱直徑包括2.54、3.81、5.08與6.35 cm。圓柱被固定在拖曳平臺(tái)上，拖曳平臺(tái)與力傳感器相連；拖曳平臺(tái)下方有2個(gè)滑輪，實(shí)驗(yàn)時(shí)，由電腦控制的拖曳平臺(tái)在搭建在水池上方的平板上做正弦式往復(fù)運(yùn)動(dòng)。只有圓柱浸入水中，其他的裝置均不會(huì)對(duì)圓柱的水動(dòng)力學(xué)結(jié)果產(chǎn)生影響。實(shí)驗(yàn)使用ATI Gamma六分力傳感器對(duì)力進(jìn)行實(shí)時(shí)記錄。本次實(shí)驗(yàn)采用直徑為3.81 cm的圓柱，與該尺寸的圓柱相比，實(shí)驗(yàn)水池范圍(1.5 m×1 m×1 m)足夠大，當(dāng)圓柱在水池中心位置處振蕩時(shí)，產(chǎn)生的渦街和尾流無(wú)法傳播到水池邊界處，因而可忽略邊界效應(yīng)，將該振蕩視為圓柱在無(wú)限大的水域中進(jìn)行的振蕩。并列圓柱的布置如圖1(c)所示，L為兩圓柱之間的間隙，D為圓柱直徑，為了便于說(shuō)明，實(shí)驗(yàn)引入間隙比(gapratio=L/D)的概念，即兩圓柱間距與圓柱直徑的比值。振動(dòng)位移隨時(shí)間的變化可表示為x=Asinωt,其中，A是振幅，ω是振動(dòng)的角速度。

1.2 實(shí)驗(yàn)數(shù)據(jù)選取

振蕩圓柱體在靜止流體中所受的力可以被分解為2個(gè)部分，第一部分為沿振動(dòng)方向的力，可以用莫里森方程表示：

第二部分為垂直于振動(dòng)方向的力(稱為升力)，可以表示為：

圖1 實(shí)驗(yàn)布置和模型Fig.1 Experiment setup and models

其中Cd、Cm和Cl分別是阻力系數(shù)，附加質(zhì)量系數(shù)和升力系數(shù)。三者主要由KC數(shù)(2πA/D)、β數(shù)(D2/Tυ)和間隙比(L/D)三個(gè)無(wú)量綱數(shù)決定，其中，T是振動(dòng)周期，υ是流體的運(yùn)動(dòng)學(xué)粘度系數(shù)。實(shí)驗(yàn)中，KC數(shù)從1～20變化，β數(shù)的選取范圍是350～2 810。雙圓柱的間隙比設(shè)定為0.5、1.0、2.0和3.0。每個(gè)實(shí)驗(yàn)工況重復(fù)3次，以保證實(shí)驗(yàn)結(jié)果的可靠性。

2 實(shí)驗(yàn)結(jié)果與分析

2.1 實(shí)驗(yàn)驗(yàn)證

并列排列，KC=12.1，β=550，L/D=0.5。

圖2給出了在靜水中，并列雙圓柱(直徑5.08 cm，間隙比0.5，KC數(shù)12.1，β數(shù)550)的位移和阻力系數(shù)Cd隨時(shí)間的變化圖中，圓柱的位移隨時(shí)間正弦變化，阻力呈現(xiàn)多頻變化，且相位落后于位移。

圖2 無(wú)量綱位移和振動(dòng)方向的 阻力系數(shù)隨著時(shí)間的關(guān)系Fig.2 non-dimensional motion and drag force vs non-dimensional time

Bearman和Sarpkaya分別對(duì)振蕩流體中的靜止圓柱進(jìn)行了類似實(shí)驗(yàn)，圖3比較了本文結(jié)果與Sarpkaya[7]的結(jié)果，可以發(fā)現(xiàn)兩組實(shí)驗(yàn)的趨勢(shì)相同，但數(shù)值上一定差異。兩組實(shí)驗(yàn)的區(qū)別在于，在Sarpkaya[7]的實(shí)驗(yàn)中，靜止圓柱被置于振蕩流場(chǎng)中，而本次實(shí)驗(yàn)是將做標(biāo)準(zhǔn)正弦振蕩的圓柱置于靜止流場(chǎng)。此外，兩組實(shí)驗(yàn)選取的β數(shù)略有差異(本文實(shí)驗(yàn)β數(shù)為1 400，Sarpkaya[7]中β數(shù)為1 380)。然而兩種實(shí)驗(yàn)方法的結(jié)果都表明，隨著KC數(shù)的增加，附加質(zhì)量系數(shù)Cm有所下降，而阻力系數(shù)Cd基本呈上升趨勢(shì)。

(本次研究β=1 400; Sarpkaya: β=1 380[7])圖3 Cd和Cm隨著KC數(shù)的變化以及 本次實(shí)驗(yàn)和Sarpkaya實(shí)驗(yàn)結(jié)果對(duì)比Fig.3 Comparison of Cdand Cm vs KC number between our experiment β=1 400 and Sarpkaya’s experiment β=1 380

2.2 阻力系數(shù)

阻力系數(shù)Cd是莫里森公式中的一個(gè)關(guān)鍵系數(shù)。與單圓柱相比，雙圓柱系統(tǒng)在流場(chǎng)中振蕩時(shí)，兩圓柱之間會(huì)產(chǎn)生相互作用，因而會(huì)使得并列雙圓柱的流場(chǎng)比單圓柱的復(fù)雜。圖4(a)和(b)分別給出了β數(shù)為550和720工況，不同間隙比情況下的平均阻力系數(shù)Cd隨著KC數(shù)的變化。隨著間隙比減小，雙圓柱平均阻力系數(shù)Cd有增大趨勢(shì)，在間隙比為0.5(試驗(yàn)中所采用的最小的間隙比)，KC數(shù)大于5時(shí)，阻力系數(shù)Cd明顯大于單圓柱振蕩時(shí)的阻力系數(shù)Cd。而在KC數(shù)較小(3,4)時(shí)，二者差異不大。這是由于正弦振蕩中的圓柱，KC數(shù)可以表示為KC=2πA/D，A是圓柱振蕩的振幅，D是圓柱的直徑，小KC數(shù)對(duì)應(yīng)的圓柱振幅較小，激發(fā)的尾流較弱，雙圓柱之間的相互作用因而也相對(duì)較弱，很難觀察到阻力系數(shù)Cd的增強(qiáng)現(xiàn)象。當(dāng)間隙比較大(3.0)時(shí)，兩者數(shù)值相近。這是由于當(dāng)兩個(gè)圓柱距離較遠(yuǎn)時(shí)，兩圓柱之間的水動(dòng)力相互作用很弱，故呈現(xiàn)出與單圓柱相一致的水動(dòng)力特性。

圖4 不同間隙比下的并列雙圓柱平均阻力系數(shù) Cd隨KC數(shù)的變化Fig.4 Drag coefficient vs KC number under different gap ratio

2.3 附加質(zhì)量系數(shù)

附加質(zhì)量系數(shù)Cm是莫里森公式中的另一個(gè)關(guān)鍵系數(shù)。由于雙圓柱振蕩引起的流場(chǎng)與單個(gè)圓柱的流場(chǎng)不一致，流體質(zhì)點(diǎn)的運(yùn)動(dòng)情況發(fā)生了明顯變化，從而導(dǎo)致附加質(zhì)量系數(shù)Cm不同。

圖5 不同間隙比下的并列雙圓柱平均 附加質(zhì)量系數(shù)隨KC數(shù)的變化Fig.5 Added mass coefficient vs KC number under different gap ratio

圖5(a)和(b)分別給出了β數(shù)等于550和720時(shí)，附加質(zhì)量系數(shù)Cm隨KC數(shù)的變化。由圖可見(jiàn)，Cm隨著KC數(shù)的增大而基本呈減小趨勢(shì)；且圓柱之間的距離越近，附加質(zhì)量系數(shù)Cm越大。與阻力系數(shù)Cd類似，當(dāng)間隙比較大(3.0)時(shí)，由于兩個(gè)圓柱相互作用弱，因此附加質(zhì)量系數(shù)Cm也與單個(gè)圓柱的情況相類似。但是在間隙比為0.5時(shí)，盡管在小KC數(shù)下兩者數(shù)值相近，但是當(dāng)KC數(shù)大于6時(shí)，0.5間隙比下雙圓柱的附加質(zhì)量系數(shù)明顯大于單圓柱的附加質(zhì)量系數(shù)，而且前者很快趨于平穩(wěn)。

同時(shí)，圖5(a)和(b)顯示，同樣的β數(shù)下，當(dāng)間隙比為1.0時(shí)，Cm約在KC數(shù)等于6時(shí)開(kāi)始趨于穩(wěn)定。當(dāng)間隙比為2.0時(shí)，Cm約在KC數(shù)等于8時(shí)趨于穩(wěn)定。當(dāng)間隙比為3.0時(shí)，Cm約在KC數(shù)等于11時(shí)趨于穩(wěn)定。單圓柱的情況下(也可以理解為間隙比為∞)，β=550時(shí)，Cm在實(shí)驗(yàn)的KC數(shù)范圍內(nèi)未出現(xiàn)穩(wěn)定，β=720時(shí)，Cm在KC數(shù)等于13時(shí)趨于穩(wěn)定。很明顯，隨著間隙比的增大，Cm在更大的KC數(shù)達(dá)到穩(wěn)定。在到達(dá)穩(wěn)定之前，Cm隨著KC的增加而減小。

2.4 升力系數(shù)

圖6為β數(shù)為550時(shí)，不同間隙比下的升力系數(shù)的均方根Cl, rms隨KC數(shù)變化的情況。

由圖可知，在KC數(shù)等于5～10的范圍內(nèi)，并列圓柱的升力系數(shù)均方根值要大于單圓柱的。同時(shí)，隨著KC數(shù)的增大，升力系數(shù)均方根值逐漸降低，但存在波動(dòng)，而且間隙比越小,該波動(dòng)對(duì)應(yīng)的KC數(shù)也越小。比如單個(gè)圓柱的波動(dòng)出現(xiàn)在KC數(shù)等于11的位置，間隙比為1.0的情況中，波動(dòng)出現(xiàn)在KC數(shù)等于7的位置，而間隙比為0.5的情況中，波動(dòng)出現(xiàn)在KC數(shù)等于6的位置。

圖6 β=550時(shí)，間隙比為0.5～4.0 下升力系數(shù)均方根隨KC數(shù)的變化Fig.6 Root mean square of lift force coefficient vs KC number (β=550, gap ratio=0.5～4.0)

為進(jìn)一步對(duì)比并列圓柱和單圓柱阻力系數(shù)Cd和升力系數(shù)Cl的不同，實(shí)驗(yàn)選取單圓柱(KC=9，β=550)，將所測(cè)升力系數(shù)Cl和阻力系數(shù)Cd的時(shí)域數(shù)據(jù)經(jīng)過(guò)傅立葉變換，得到Cl和Cd的頻譜(見(jiàn)圖7)?？梢钥吹?，兩個(gè)系數(shù)的一階頻均在0.77處，同時(shí)升力系數(shù)Cl在1.55(約二倍頻率處)處出現(xiàn)了明顯的二階分量，且二階分量占優(yōu)。阻力系數(shù)Cd在2.33處出現(xiàn)較明顯的二階分量(約三倍頻率處)，但一階分量占優(yōu)。

此外，對(duì)β=550條件下，不同KC數(shù)(1～20)的單圓柱和雙圓柱(L/D=0.5)的實(shí)驗(yàn)數(shù)據(jù)進(jìn)行傅立葉變換，得到x坐標(biāo)為階頻數(shù)，y坐標(biāo)為KC數(shù)，z坐標(biāo)為譜密度的3D頻譜圖(見(jiàn)圖8、9)。

圖7 單圓柱(KC=9,β=550) 情況下Cl和Cd的傅里葉頻譜Fig.7 Fourier spectrum for single cylinder (KC=9, β=550)

圖8 單個(gè)圓柱情況下傅里葉變化 得到的Cl頻譜圖

圖9 雙圓柱并列振蕩情況下(間隙比0.5)傅里葉 變化得到的Cl頻譜圖Fig.9 Fourier spectrum for dual cylinders’ lift coefficient (gap ratio=0.5)

對(duì)比圖8、9，可以看到單圓柱情況下，隨著KC數(shù)的增大，圓柱的升力系數(shù)頻率先是一階占優(yōu)，而后逐漸過(guò)渡到二階占優(yōu)(KC=11)。并列圓柱情況下，升力系數(shù)主要表現(xiàn)為二階頻率占優(yōu)，并且存在較為明顯的四階頻率分量。

3 數(shù)值模擬結(jié)果與分析

由于在水槽中進(jìn)行的實(shí)驗(yàn)很難將流場(chǎng)清晰地展現(xiàn)出來(lái)，因此，為了進(jìn)一步分析雙圓柱并列振蕩和單圓柱振蕩的流場(chǎng)特征，本次研究開(kāi)展了基于BDIM(Boundary Data Immersion Method)的數(shù)值模擬，得到不同振蕩形式下的流場(chǎng)實(shí)時(shí)數(shù)據(jù)和流場(chǎng)的渦量云圖。研究將通過(guò)數(shù)值模擬實(shí)時(shí)數(shù)據(jù)計(jì)算升力系數(shù)、阻力系數(shù)和附加質(zhì)量系數(shù)，并設(shè)置控制體，從能量的角度解釋水槽實(shí)驗(yàn)中阻力系數(shù)增大的現(xiàn)象。Zhao[10]的研究證明，與三維模型相比，二維數(shù)值模型也能夠清晰的描繪出流場(chǎng)特征并準(zhǔn)確的計(jì)算圓柱的水動(dòng)力學(xué)載荷，因此，本次研究使用二維模型進(jìn)行模擬。

數(shù)值模擬的計(jì)算中將所有物理量都無(wú)量綱化：(x,y)=(x’,y’)/D、 (u,v)=(u’,v’)/Um、p=p’/ρUm、t=Umt’/D。 (x,y)是網(wǎng)格上的坐標(biāo)。 (u,v)是x、y方向上的速度分量t、p、ρ和Um分別是時(shí)間，壓強(qiáng)，流體密度和圓柱的最大振蕩速度，則無(wú)量綱化的控制方程化為：

數(shù)值模擬布置如圖10所示，虛線表示控制體的邊界?？刂企w長(zhǎng)度為8D，寬度為4D，以滿足大運(yùn)動(dòng)幅度的要求。因?yàn)閿?shù)值模擬的邊界是不可穿透的，因此將計(jì)算域的大小設(shè)置為40D×20D，以達(dá)到消除邊界效應(yīng)和阻擋效應(yīng)。 Lily-Pad的求解器不需要我們建立網(wǎng)格，屏幕的像素將直接成為網(wǎng)格的節(jié)點(diǎn)并建立非結(jié)構(gòu)化網(wǎng)格。然而，求解器要求必須以像素為單位定義柱面和計(jì)算域的大小，因此在計(jì)算中，為了平衡計(jì)算效率和精度，每個(gè)圓柱體的直徑D被設(shè)置為50個(gè)像素。

圖10 數(shù)值模擬布置圖Fig.10 Numerical simulation setup

4.1 正確性和收斂性驗(yàn)證

首先采用三種不同全局精度的網(wǎng)格(delta/D=0.025、0.02和0.016)對(duì)振蕩單圓柱繞流進(jìn)行數(shù)值模擬，以檢驗(yàn)網(wǎng)格的收斂性。收斂性驗(yàn)證中，KC數(shù)目設(shè)置為6，雷諾數(shù)Re設(shè)置為120。計(jì)算所得的附加質(zhì)量系數(shù)和阻力系數(shù)見(jiàn)表1。以delta/D=0.016所得到的Cm值(1.975)和Cd值(1.979)為參考，計(jì)算不同網(wǎng)格精度條件下，附加質(zhì)量系數(shù)Cm和阻力系數(shù)Cd的相對(duì)誤差和。

由表1可知，與delta/D=0.016的工況相比，delta/D=0.02工況的附加質(zhì)量系數(shù)的相對(duì)誤差為0.2%，阻力系數(shù)的相對(duì)誤差為0.56%。相比之下，delta/D=0.025工況的附加質(zhì)量系數(shù)的相對(duì)誤差為6.26%，阻力系數(shù)的相對(duì)誤差為8.38%，遠(yuǎn)大于delta/D=0.02工況的相對(duì)誤差，證實(shí)了delta/D=0.02的網(wǎng)格滿足收斂性要求。圖11(b)和(c)進(jìn)一步給出了delta/D為0.02和0.016時(shí)，圓柱向右以最大速度通過(guò)平衡位置的的渦量云圖(其相位定義為π/2)，可以看到兩者幾乎相同。然而，圖11(a)中delta/D為0.025的結(jié)果明顯不同，流線對(duì)稱，未觀察到圓柱兩側(cè)交替泄渦現(xiàn)象。

為了證明數(shù)值模擬算法的正確性，在進(jìn)行數(shù)值模擬研究前，研究首先使用前人在精度更高的網(wǎng)格上計(jì)算得到的結(jié)果進(jìn)行比對(duì)：另取與Tong[9]相同的工況參數(shù)設(shè)置，即KC數(shù)為5、雷諾數(shù)Re為100，并將該工況下使用三種精度更高的網(wǎng)格計(jì)算得到的Cm和Cd的結(jié)果和Tong[9]進(jìn)行對(duì)比，見(jiàn)表2。

表1 網(wǎng)格收斂性驗(yàn)證Table 1 Grid convergence verification

表2 Tong(2015)[9]數(shù)值模擬參數(shù)和本次數(shù)模擬參數(shù)的對(duì)比Table 2 Comparison between Tong (2015)[9] numerical simulation parameters and current simulation parameters

Tong[9]設(shè)置的三種網(wǎng)格delta/D分別為0.003, 0.001, 0.000 75，其精度較本次研究使用的網(wǎng)格更高，其計(jì)算得到的(Cd，Cm)分別為(2.11，2.43)、(2.11，2.43)、(2.11，2.42)、最終Tong[9]選取了網(wǎng)格2即delta/D為0.001的網(wǎng)格進(jìn)行數(shù)值模擬。與網(wǎng)格2相比對(duì)，本次研究中網(wǎng)格4計(jì)算的結(jié)果(2.29，2.23)與之差距較大，網(wǎng)格5得到的Cd和Cm略小于網(wǎng)格2計(jì)算的結(jié)果，網(wǎng)格6計(jì)算的Cd小于網(wǎng)格2計(jì)算的Cd，而Cm略大于網(wǎng)格2。同時(shí)，與網(wǎng)格6相比，網(wǎng)格5計(jì)算的Cd更接近于Tong[9]使用的精度更高的網(wǎng)格2的Cd值。因此網(wǎng)格5和網(wǎng)格6滿足收斂性驗(yàn)證，且其計(jì)算結(jié)果與前人在精度更高的網(wǎng)格上的計(jì)算結(jié)果更為接近，滿足正確性驗(yàn)證。同時(shí)，從網(wǎng)格5和網(wǎng)格6中都可以觀察到交替泄渦的現(xiàn)象。為了提高運(yùn)算效率，研究采取網(wǎng)格5(delta/D=0.02)進(jìn)行數(shù)值模擬計(jì)算。

3.2 阻力系數(shù)Cd分析

在并列雙圓柱的數(shù)值模擬中，雷諾數(shù)Re選定為120，KC數(shù)為6，β數(shù)(Re/KC)為20，間隙比選定為0.50、0.75、1.00、1.25、1.50、2.00和3.00，同時(shí)選取單個(gè)圓柱振蕩進(jìn)行對(duì)比(見(jiàn)圖12)，圖中的虛線所對(duì)應(yīng)的縱坐標(biāo)值代表同一KC數(shù)和β數(shù)下的單個(gè)圓柱阻力系數(shù)Cd的值。

圖11 在KC=6時(shí)間為65時(shí)，不同網(wǎng)格的渦旋場(chǎng)的比較Fig.11 Comparison of vortices in different grids when KC=6, non-dimensional time=65

圖12 KC=6時(shí)的阻力系數(shù)CdFig.12 Drag Coefficient (KC=6)

可見(jiàn)，振蕩并列圓柱的阻力系數(shù)Cd比單圓柱的較大，隨著間隙比的增大，阻力系數(shù)Cd單調(diào)減小，并趨近于單圓柱的阻力系數(shù)。這與實(shí)驗(yàn)結(jié)果的趨勢(shì)一致。

4 數(shù)值模擬的流場(chǎng)分析

4.1 振蕩雙圓柱流場(chǎng)

圖13給出了間隙比為0.5振蕩雙圓柱繞流在一個(gè)周期內(nèi)的4個(gè)特征時(shí)刻的流場(chǎng)，在這4個(gè)時(shí)刻圓柱的運(yùn)動(dòng)狀態(tài)分別表示左側(cè)最大位移處(定義相位為0)，向右通過(guò)平衡位置處(定義相位為π/2)，右側(cè)最大位移處(定義相位為π)，向右通過(guò)平衡位置處(定義相位為3π/2)，以此描述一個(gè)運(yùn)動(dòng)周期內(nèi)的流場(chǎng)變化。圖13(a)是圓柱運(yùn)動(dòng)到最左側(cè)位置(相位為0)時(shí)的流場(chǎng)，此時(shí)圓柱泄出渦對(duì)A和B，A的下部渦流(藍(lán)色的逆時(shí)針?shù)鰷u)和B的上部渦流(紅色的順時(shí)針?shù)鰷u)攜帶了大量能量。并在間隙中產(chǎn)生強(qiáng)相互作用，形成強(qiáng)噴射流X，快速遠(yuǎn)離圓柱，其速度超過(guò)了圓柱的振蕩速度。當(dāng)圓柱移動(dòng)到右側(cè)最大位移時(shí)，如圖13(c)所示，渦對(duì)C和D泄出，與A和B類似，C的下部渦流和D的上部渦流在間隙中形成了新的渦流對(duì)Y。與X類似，Y攜帶了大量能量并快速遠(yuǎn)離圓柱。同時(shí)圓柱激發(fā)了新的渦流對(duì)E與F。隨后，如圖13(d)所示，此后圓柱將運(yùn)動(dòng)到左側(cè)最大位移處，間隙中產(chǎn)生新的射流渦對(duì)。這些攜帶能量且具有較快速度的渦對(duì)，是雙圓柱系統(tǒng)阻力系數(shù)增大的原因。

圖13 雙圓柱的流場(chǎng)發(fā)展 (KC數(shù)為6, 間隙比為0.5)

5.2 平均能量流出率

為了從能量角度定量的解釋雙圓柱系統(tǒng)阻力系數(shù)增大的現(xiàn)象，本文引入平均能量流出率(Mean Energy Outflow Rate，簡(jiǎn)稱)，其定義式為：

其中：u為無(wú)量綱流體速度，其方向是控制體的外法線方向；ds是控制體的無(wú)量綱的單位邊長(zhǎng)，控制體設(shè)置示意圖見(jiàn)圖10。平均能量流出率的物理含義是指：一個(gè)完整的振蕩周期中，選定控制體在單位時(shí)間內(nèi)能量的流出，其大小體現(xiàn)了該情況下由于圓柱振蕩而流出的能量多少。選取KC數(shù)等于6的情況，對(duì)不同間隙比下的雙圓柱系統(tǒng)求得平均能量流出率(控制體內(nèi)包含圓柱體)，如圖14所示，其中水平虛線表示振蕩單圓柱工況的平均能量流出率。

圖14 KC數(shù)等于6時(shí)不同 間隙比下圓柱求得平均能量流出率Fig.14 Energy outflow rate of cylinder with different gap ratio(KC=6)

因?yàn)槠骄芰苛鞒雎适菍?duì)圓柱相對(duì)于水的相對(duì)速度進(jìn)行邊界積分，而相對(duì)速度的大小影響了莫里森公式中的阻力，系數(shù)Cd，因此平均能量流出率隨間隙比的變化從能量的角度解釋了小間隙比條件下阻力系數(shù)Cd增大的機(jī)制。

5 結(jié)語(yǔ)

本文對(duì)做正弦式往復(fù)振蕩的并列雙圓柱系統(tǒng)在靜止流場(chǎng)里的水動(dòng)力性能進(jìn)行了實(shí)驗(yàn)和數(shù)值模擬。實(shí)驗(yàn)發(fā)現(xiàn)，小間隙比并列振蕩的情況與單個(gè)圓柱相比，阻力系數(shù)Cd在KC數(shù)>5時(shí)出現(xiàn)明顯的增大現(xiàn)象，附加質(zhì)量系數(shù)Cm在小KC數(shù)下提前趨于穩(wěn)定，升力系數(shù)Cl的方均根增大，而且Cl在小KC數(shù)下提前出現(xiàn)二階頻率分量。隨后，使用基于BDIM的數(shù)值模擬，針對(duì)KC數(shù)為6，β數(shù)為20(Re=120)的情況進(jìn)行水動(dòng)力學(xué)載荷計(jì)算和流場(chǎng)分析。數(shù)值模擬通過(guò)使用delta/D=50的網(wǎng)格，發(fā)現(xiàn)雙圓柱系統(tǒng)在小間隙比下并列振蕩時(shí)，由于兩圓柱流場(chǎng)之間的相互作用，在間隙中會(huì)產(chǎn)生間隙流和攜帶能量的渦流。通過(guò)計(jì)算各系統(tǒng)中流出的能量，研究發(fā)現(xiàn)小間隙比下雙圓柱系統(tǒng)產(chǎn)生的渦流攜帶的能量高于單圓柱系統(tǒng)，這一攜帶更多能量的渦流導(dǎo)致了小間隙比下雙圓柱系統(tǒng)阻力系數(shù)Cd的增大，而較大間隙比下雙圓柱之間的作用比較微弱，流出的能量也接近與單圓柱振蕩流出的能量。