小間距比下串列雙圓柱渦激振動數(shù)值模擬研究：振動響應(yīng)和流體力

陳威霖， 及春寧， 許 棟

(天津大學(xué) 水利工程仿真與安全國家重點實驗室，天津 300072)

圓柱渦激振動是流固耦合的經(jīng)典問題，近幾十年對該現(xiàn)象的研究取得了大量的成果[1-3]。此外，渦激振動也常見于海洋、航空航天以及核能等眾多的工程領(lǐng)域，比如輸油立管、冷卻管等，對該現(xiàn)象研究具有非常重要的實際意義。影響渦激振動的眾多因素中，雷諾數(shù)和質(zhì)量比是最受關(guān)注的。在低雷諾數(shù)(～O(102))時，圓柱的響應(yīng)呈現(xiàn)出初始和下端兩個分支[4-5],其中最大振幅出現(xiàn)在下端分支上；此時，質(zhì)量比和阻尼不再是影響圓柱響應(yīng)的主要因素。當雷諾數(shù)較高(≥O(103))時，質(zhì)量-阻尼系數(shù)成為衡量圓柱響應(yīng)的一個重要因素[6-7]：當質(zhì)量-阻尼系數(shù)較低時，響應(yīng)呈現(xiàn)出三個分支，分別為初始、上端和下端分支，其中最大振幅出現(xiàn)在上端分支上；當質(zhì)量-阻尼系數(shù)較高時，響應(yīng)則僅為兩個分支，分別為初始和下端分支，其中最大振幅出現(xiàn)在初始分支上。

當圓柱的數(shù)量不再限定為一個，圓柱與圓柱、圓柱與流體之間的耦合作用變得更加復(fù)雜[8-9]。串列布置為眾多情形中比較常見的一種，圓柱之間的相對位置決定了下游圓柱始終要受到上游圓柱脫落旋渦的影響，而下游圓柱對上游圓柱的影響則只存在于一定的間距范圍內(nèi)[10]。Prasanth等[11]對低雷諾數(shù)Re=100和大間距比L*=L/D=5.5(其中L為兩圓柱中心的距離)下串列雙圓柱渦激振動進行了數(shù)值模擬研究，發(fā)現(xiàn)下游圓柱響應(yīng)要明顯大于對應(yīng)的單圓柱渦激振動，但類似于高雷諾數(shù)下單圓柱渦激振動的情況。受到上游圓柱的影響，下游圓柱的鎖定區(qū)間要明顯大于單圓柱渦激振動情況；此外，研究還發(fā)現(xiàn)大振幅的渦激振動不僅出現(xiàn)在鎖定區(qū)域，在非鎖定區(qū)間，響應(yīng)振幅也可以很大。Zhao等[12]對直徑比為0.5的串列雙圓柱渦激振動在間距比L*=1.5～3.0范圍內(nèi)展開了研究，其中雷諾數(shù)為Re=200。研究發(fā)現(xiàn)，當間距比為1.5～2.0時，圓柱響應(yīng)出現(xiàn)四個不同的分區(qū)。Bao等[13]應(yīng)用基于有限元的數(shù)值方法研究了不同頻率比下串列雙圓柱渦激振動的情況，其中間距比L*=5.0和Re=150。研究表明，此時下游圓柱對上游圓柱的影響幾乎可以忽略，而上游圓柱對下游圓柱則有顯著的影響。此外，下游圓柱的流向響應(yīng)對頻率比的敏感程度要明顯高于橫向的情況。Mysa等[14]研究了串列雙圓柱渦激振動中對響應(yīng)起關(guān)鍵作用的因素，其中間距比L*=4.0～10.0和雷諾數(shù)Re=100。研究發(fā)現(xiàn)，上游圓柱的尾流與下游圓柱發(fā)生相互作用，其中與圓柱運動同相的成分對下游圓柱的振動起到了促進作用。King等[15]通過水槽實驗研究了Re=103～2×104下串列雙圓柱渦激振動的情況，其中間距比為L*=1.25～7.0。 結(jié)果顯示，當間距比為6.5時，上游圓柱的響應(yīng)與對應(yīng)雷諾數(shù)下的單圓柱渦激振動相似。受上游圓柱的影響，下游圓柱的大振幅響應(yīng)一直保持到實驗的最大折合流速。Brika等[16]對雷諾數(shù)Re=5×103～2.7×104下串列雙圓柱渦激振動進行了一系列的研究，其中間距比為L*=7～25。隨著間距比的增加，上游圓柱的響應(yīng)越來越接近單圓柱渦激振動的情況；即使在間距比L*=16～25下，上游圓柱對下游圓柱的影響仍然存在。Assi等[17]通過水槽實驗研究了上游圓柱固定情況下的串列雙圓柱渦激振動情況，相應(yīng)的雷諾數(shù)為Re=3×103～1.3×104，間距比為L*=2.0～5.6. 在間距比L*=3.0～5.6時，下游圓柱的振幅隨著折合流速增加而持續(xù)增加，出現(xiàn)了尾流弛振現(xiàn)象。此后，Assi等[18]對該尾流弛振現(xiàn)象通過一系列的實驗進行更深入的研究，發(fā)現(xiàn)尾流弛振來源于下游圓柱與上游圓柱尾流之間的不穩(wěn)定的旋渦-結(jié)構(gòu)之間的耦合作用。

目前對串列雙圓柱渦激振動的研究已經(jīng)取得了相當多的成果，但是對小間距比(L*=1.1～1.5)下該現(xiàn)象的研究仍較少。且本文模擬結(jié)果顯示，在小間距比下，圓柱之間存在較強的耦合作用，而且發(fā)現(xiàn)了一些新的現(xiàn)象，比如低雷諾數(shù)下串列雙圓柱的尾流弛振現(xiàn)象等。

1 數(shù)值方法

1.1 控制方程

流體運動的數(shù)值模擬采用浸入邊界法[19]，無量綱的控制方程如下

(1)

(2)

對以上控制方程采用二階精度的Adams-Bashforth時間格式進行離散，可得控制方程的守恒形式如下

(3)

(4)

針對傳統(tǒng)浸入邊界法施加邊界條件精度不高的情況，Ji等[19]提出了基于嵌入式迭代的浸入邊界法，將浸入邊界法嵌入到壓強泊松方程的迭代求解中，利用壓強的中間解比初始值更接近真實值的特點，迭代修正附加體積力，在不顯著增加計算耗時的前提下，提高整個算法的求解精度。

對僅做橫流向運動的剛性圓柱體，其無量綱運動方程可以用下述方程來描述

(5)

式中:y為無量綱位移;t為無量綱時間;m*為圓柱體質(zhì)量比;FN=fnD/U∞為無量綱固有頻率(fn為固有頻率);ζ為結(jié)構(gòu)阻尼比;Cl為圓柱受到的橫流向流體力系數(shù)。方程采用標準的Newmark-β法求解。

以上流體和剛體運動控制方程采用圓柱直徑D、來流速度U∞和流體密度ρ進行無量綱化。

1.2 問題描述

擬對小間距比下串列雙圓柱渦激振動展開深入研究，數(shù)值模擬相關(guān)參數(shù)設(shè)置如下：兩圓柱的間距比L*=L/D=1.1～1.5，Re=U∞D(zhuǎn)/υ=100(其中υ為運動黏性系數(shù));折合流速Ur=U∞/fnD=3～30(最小間隔為ΔUr=0.5);質(zhì)量比m*=4m/ρπD2=2.0(其中m圓柱的質(zhì)量)。

為消除邊界對圓柱渦激振動的影響，本文采用了較大的計算域,如圖1所示。將計算域的坐標原點設(shè)在串列兩圓柱距離的中心上，距離入口和出口邊界為100D，距離上和下邊界為50D。因此，相應(yīng)的阻流比為B=D/H=0.01。Sen等[20]對圓柱繞流數(shù)值模擬研究發(fā)現(xiàn)，當阻流比B≤0.01時，阻流比對流動特征參數(shù)的影響不再重要。

圖1 計算域與邊界條件

采用正交笛卡爾網(wǎng)格對計算區(qū)域劃分，X(流向)和Y(橫向)方向的網(wǎng)格數(shù)分別為1 024和512，其中加密區(qū)域內(nèi)無量綱網(wǎng)格尺寸Δx/D=Δy/D=1/64。加密區(qū)域為6.5D×4.0D;編號1和2分別為上游和下游圓柱。

計算域邊界條件設(shè)置如下。入口邊界為Dirichlet型邊界條件(u=U∞,ν=0)；出口邊界為Neumann型邊界條件(?u/?x=0,?v/?x=0)；上下邊界為自由滑移邊界條件。

為了滿足Courant-Friedrichs-Lewy條件，CFL數(shù)滿足CFL=UmaxΔt/Δx≤0.5，其中，Umax為最大流速，Δt為時間步長,取為Δt=0.006。

2 方法驗證

本文采用的數(shù)值模型和程序已從單圓柱繞流[21-22]，單圓柱渦激振動以及串列并列雙圓柱[23-25]等多角度進行了驗證，讀者可自行查閱。

3 結(jié)果和討論

3.1 振動響應(yīng)

當間距比為L*=1.1，兩圓柱的響應(yīng)從折合流速Ur=3.0開始增加，當折合流速Ur>4.0以后，圓柱的響應(yīng)開始急劇增加，并在折合流速Ur=9.0時兩圓柱響應(yīng)同時達到最大值，分別為Yrms/D=0.39和Yrms/D=0.58。上游圓柱的最大響應(yīng)與單圓柱渦激振動的最大振幅[27]Yrms/D=0.40非常接近，而下游圓柱則要明顯大于單圓柱渦激振動的情況。但是，相比于較大間距比下的串列雙圓柱渦激振動[28-29]，如圖2(a)所示，兩圓柱的響應(yīng)則小了一些。兩圓柱的響應(yīng)均在Ur=9.0以后發(fā)生下降，上游圓柱分別在折合流速Ur=12.0和Ur=14.0時經(jīng)歷一個谷值和一個峰值，其大小為Yrms/D=0.27和Yrms/D=0.39；下游圓柱同樣經(jīng)歷了一個谷值(Ur=12.0)和一個峰值(Ur=20.0)，其大小為Yrms/D=0.35和Yrms/D=0.60。此后，下游圓柱的響應(yīng)隨折合流速急劇下降，而上游圓柱則經(jīng)歷了一個比較平緩的區(qū)域(Ur=17.0～22.0)才急劇下降。

從圓柱響應(yīng)隨折合流速變化的情況可知，此時圓柱響應(yīng)不再與單圓柱渦激振動相似，沒有出現(xiàn)初始和下端分支[30]，而且圓柱的大振幅響應(yīng)出現(xiàn)在更大的折合流速范圍內(nèi)(Ur=4.0～28.0)，遠大于單圓柱渦激振動的范圍(Ur=4.0～8.5)。這些不同反映了串列雙圓柱渦激振動在小間距比下強烈的耦合作用，而該耦合機制將在下部分進行詳細的分析。

需要說明的是在折合流速從Ur=3.5～4.0的過程中出現(xiàn)了一個突然的下降，這是圓柱與流體之間相互作用的結(jié)果，此現(xiàn)象的機制已在文獻[27]中提及。

當間距比為L*=1.2時，以折合流速Ur=4.0和Ur=15.0為界，串列雙圓柱的響應(yīng)呈現(xiàn)出三個分支。圓柱響應(yīng)在Ur=3.0～4.0時緩慢增加直到Ur>4.0以后響應(yīng)突然增加，圓柱響應(yīng)振幅分別在折合流速Ur=14.0和Ur=12.0時取得最大值Yrms/D=0.44和Yrms/D=0.63，該最大值要稍大于間距比L*=1.1的情況。在折合流速Ur>15.0以后，兩圓柱的振幅出現(xiàn)急劇下降，并均在折合流速Ur=20.0時達到一個谷值(Yrms/D=0.29)；此后，下游圓柱響應(yīng)隨著折合流速而持續(xù)增加，而上游圓柱則要經(jīng)歷一個平穩(wěn)區(qū)域(Ur=20.0～28.0)，之后圓柱響應(yīng)隨折合流速增加。在大折合流速下，圓柱響應(yīng)隨折合流速增加的情況類似于高雷諾數(shù)下串列雙圓柱渦激振動的尾流弛振現(xiàn)象[31-33]，為首次出現(xiàn)在低雷諾數(shù)串列雙圓柱渦激振動中；說明了該間距比下圓柱與流體之間的耦合作用與此前不再相同，其耦合機制將在下部分詳細闡述。

(a) L/D=1.1

(b) L/D=1.2

(c) L/D=1.3

(d) L/D=1.5

當間距比為L*=1.3時，串列雙圓柱的響應(yīng)均呈現(xiàn)出四個分支。第一分支出現(xiàn)在Ur=3.0～4.0范圍內(nèi)，此時圓柱響應(yīng)稍大于零且不隨折合流速變化；當折合流速Ur>4.0以后，圓柱響應(yīng)進入第二分支，此時振幅急劇增加，一直到Ur=8.5以后，圓柱響應(yīng)出現(xiàn)另一個更急劇的上升(第三分支)，兩圓柱響應(yīng)最大值Yrms/D=0.46和Yrms/D=0.74同時在Ur=9.0時取得。該最大值要稍大于更小間距比L*=1.1～1.2的情況。在折合流速Ur=9.0～23.0內(nèi)，上游圓柱響應(yīng)幾乎保持不變而下游圓柱響應(yīng)隨折合流速緩慢下降。之后，兩圓柱響應(yīng)均急劇下降，在折合流速Ur=24.0時達到谷值(Yrms/D=0.17)；圓柱響應(yīng)此后隨折合流速增加而持續(xù)增加(第四分支)。需說明的是此時下游圓柱的響應(yīng)要明顯小于上游圓柱的情況，該現(xiàn)象僅在此間距比下出現(xiàn)，此時圓柱與流體之間的耦合作用也與其余情況不同。

當間距比為L*=1.5時，圓柱響應(yīng)與間距比L*=1.3的情況相同在折合流速Ur=3.0～4.0內(nèi)在零附近且保持不變；之后圓柱響應(yīng)隨折合流速快速增加，上游圓柱在折合流速Ur=6.5時取得最大值Yrms/D=0.54，下游圓柱的最大值Yrms/D=0.77則在Ur=7.5時得到。之后，圓柱響應(yīng)隨折合流速增加而緩慢減小，在折合流速Ur≥28.0之后，圓柱響應(yīng)不再變化。在模擬的最大折合流速Ur=30.0下，兩圓柱的響應(yīng)振幅分別為Yrms/D=0.35和Yrms/D=0.55，均要遠大于單圓柱渦激振動的情況，而且此時下游圓柱的振幅甚至要比單圓柱渦激振動的最大振幅(Yrms/D=0.40)大37.5%。

需要說明的是上游圓柱響應(yīng)在Ur=4.5～5.0之間出現(xiàn)了類似間距比L*=1.1時的突然下降情況，但與之不同的是此時下游圓柱的響應(yīng)并未出現(xiàn)下降，該現(xiàn)象的原因分析如下。如圖3所示，當折合流速Ur=4.5時，上游圓柱和下游圓柱位移的相位差維持在φy=101°上；當折合流速Ur=5.0時，上游圓柱和下游圓柱位移的相位差維持在φy=77.5°上。對于前者，當上游圓柱位于平衡位置時，下游圓柱恰好位于向最大負位移運動的過程中，此時上游圓柱后有充足的空間使得旋渦向下游移動，因此，為上游圓柱的運動提供了足夠的動力；對于后者，當上游圓柱位于平衡位置時，下游圓柱恰好位于由最大負位移向平衡位置運動的過程中，下游圓柱的運動使得上游圓柱后的旋渦自由脫落的空間越來越小，因此，在一定程度上阻礙了上游圓柱的振動，使得其振幅下降。

圖3 折合流速Ur=4.5和Ur=5.0下串列雙圓柱位移相位差

Fig.3 Phase difference between the displacements of two tandem cylinders atUr=4.5 andUr=5.0

總結(jié)起來，如圖2所示，小間距比下串列雙圓柱渦激振動與大間距比(L*=5.0)下完全不同；此時串列雙圓柱的鎖定區(qū)間要明顯大于大間距比下的情況。在折合流速較小時(Ur≈10～15)，小間距比下串列雙圓柱的響應(yīng)呈現(xiàn)出隨折合流速增加而增加；當折合流速較大時，各間距比下串列雙圓柱的響應(yīng)則明顯不同。當L*=1.1時，串列雙圓柱的響應(yīng)呈現(xiàn)出先增加再急劇下降到接近于零，而L*=1.2時，兩圓柱的響應(yīng)隨折合流速增加而增加，并未出現(xiàn)下降。與前兩個間距比不同的是間距比L*=1.3的不規(guī)律振動出現(xiàn)在更大的折合流速下，響應(yīng)隨折合流速的增加而增加，未出現(xiàn)下降。當L*=1.5時，串列雙圓柱的響應(yīng)在大折合流速以后隨著折合流速緩慢下降。

3.2 流體力系數(shù)

(a) L/D=1.1

(b) L/D=1.2

(c) L/D=1.3

(d) L/D=1.5

Fig.4 The mean drag coefficients versus the reduced velocity at different spacing ratios

(a) L/D=1.1

(b) L/D=1.2

(c) L/D=1.3

(d) L/D=1.5

Fig.5 The RMS lift coefficients versus the reduced velocity at different spacing ratios

如圖5所示，各間距比下，串列雙圓柱的升力均方根隨折合流速呈現(xiàn)出不規(guī)律的變化。當間距比L*=1.1時，上游和下游圓柱的升力均方根呈現(xiàn)出相似的變化。上游圓柱呈現(xiàn)出四個先增后減的趨勢，折合流速Ur=3.5、Ur=4.5、Ur=6.5、Ur=15.0時分別對應(yīng)相應(yīng)的峰值，而折合流速Ur=4.0、Ur=6.0、Ur=7.0時分別對應(yīng)相應(yīng)的谷值。上游圓柱的升力均方根在Ur≥15.0以后，迅速下降并在Ur≥28.0以后穩(wěn)定下來。下游圓柱的升力均方根呈現(xiàn)出三個先增后減的趨勢，折合流速Ur=3.5、Ur=4.5、Ur=10.0時分別對應(yīng)相應(yīng)峰值，而折合流速Ur=4.0、Ur=7.0時分別對應(yīng)相應(yīng)谷值。此后下游圓柱的升力均方根隨折合流速下降，但在折合流速Ur≥22.0以后，升力均方根隨折合流速緩慢增加。

當間距比L*=1.3時，兩圓柱的升力均方根趨勢與間距比L*=1.2時的情況幾乎一致，僅在折合流速Ur≥24.0以后不再相同。在折合流速Ur=24.0～30.0內(nèi)，兩圓柱的升力均方根均較小，且隨折合流速僅稍微增加。需要說明的是，雖然此區(qū)域內(nèi)上游圓柱的響應(yīng)要大于下游圓柱，但上游圓柱的升力均方根卻要小于下游圓柱的。

值得注意的是在串列雙圓柱渦激振動中，兩圓柱的升力均方根均在折合流速Ur=7.0～8.0范圍內(nèi)出現(xiàn)一個谷值，而且在間距比L*≤1.3時，兩圓柱的升力均方根值甚至要小于振幅幾乎為零的折合流速Ur=3.0時的情況。出現(xiàn)這種現(xiàn)象的原因是該谷值對應(yīng)的折合流速下，兩圓柱的升力不再由單一頻率成分組成，其中上游圓柱的高倍頻占主導(dǎo)而下游圓柱的基頻占主導(dǎo)。以間距比L*=1.2為例，對兩圓柱升力(Ur=7.0)的小波分析表明，如圖6(a)所示，上游圓柱的升力主頻為3倍頻(f=0.429)；而下游圓柱的主頻為基頻(f=0.143)，如圖6(b)所示。此外，提取與兩圓柱運動同相的升力成分發(fā)現(xiàn)，此時與兩圓柱運動同相的升力成分較大；這說明高頻成分的存在提高了與兩圓柱運動同相的升力成分，因此保證了圓柱能以較大的振幅振動。

(a) 上游圓柱

(b) 下游圓柱

圖6 對應(yīng)間距比L*=1.2和折合流速Ur=7.0時上游和下游圓柱的升力和頻域的歷時曲線

Fig.6 Time history of the lift coefficient and the frequency of upstream and downstream cylinder atL*=1.2 andUr=7.0

在圓柱渦激振動中，通過快速傅里葉變換(FFT)對圓柱的升力脈動值進行分析得到其主導(dǎo)頻率并將其定義為St數(shù)。當間距比L*≤1.3時，由于圓柱之間的耦合作用較強，使得某些折合流速下圓柱升力的頻譜混亂，無法得到相應(yīng)的St數(shù)，比如L*=1.1的Ur=12.0～25.0；因此，本文僅考慮間距比L*=1.5時的情況。如圖7所示，兩圓柱的St數(shù)在絕大多數(shù)折合流速下是相等的。當折合流速Ur=5.5時，兩圓柱的St數(shù)不再相等，且兩頻率不是倍數(shù)關(guān)系；這反映出兩圓柱之間的耦合作用無法實現(xiàn)穩(wěn)定。該折合流速下，兩圓柱的響應(yīng)恰好處在由小振幅向大振幅的過渡，因此，可以看作是一個調(diào)制的階段。

圖7 間距比L*=1.5時St數(shù)隨折合流速變化情況

另外，在折合流速Ur=7.5～8.5范圍內(nèi)，上游和下游圓柱的主導(dǎo)頻率為基頻的三倍，而該現(xiàn)象的出現(xiàn)與圓柱升力與位移相位差的跳躍有關(guān)。

4 結(jié) 論

本文對小間距比下串列雙圓柱渦激振動進行了廣參數(shù)空間的數(shù)值模擬研究，其中Re=100，間距比為L*=1.1～1.5，折合流速為Ur=3.0～30.0。為保持兩圓柱間距不變，兩圓柱僅作橫向振動。對圓柱響應(yīng)的研究發(fā)現(xiàn)，在小間距比下，圓柱的響應(yīng)呈現(xiàn)出完全不同于大間距比下的情況；上游圓柱的響應(yīng)也與單圓柱渦激振動有很大區(qū)別，沒有初始和下端分支。根據(jù)響應(yīng)的不同，可以分為三種。當間距比L*≤1.1時，響應(yīng)存在于較大的折合流速范圍內(nèi)(Ur=4.0～28.0)，且在大折合流速內(nèi)，兩圓柱之間的耦合作用不再穩(wěn)定。當間距比L*=1.2～1.3時，響應(yīng)在小折合流速時呈現(xiàn)出類似于單圓柱渦激振動的現(xiàn)象，而類似高雷諾數(shù)串列雙圓柱中的尾流弛振現(xiàn)象出現(xiàn)在大折合流速下，且隨間距比的增加，開始該現(xiàn)象的折合流速增大。當間距比L*≥1.5時，響應(yīng)隨折合流速增加到最大值以后緩慢減小，并最終穩(wěn)定在某個值上。

在小間距比下，串列雙圓柱渦激振動的響應(yīng)幅值相比于大間距比下要小一些。上游和下游圓柱的最大響應(yīng)均隨間距比增加而增加；當間距比L*=1.5時，下游圓柱的最大振幅達到了Yrms/D=0.77，且隨間距比的進一步增加而下游圓柱的最大振幅減小，說明上游和下游圓柱之間的耦合作用在間距比L*=1.5時達到最大。

串列雙圓柱的脫渦頻率(St數(shù))反映的是圓柱之間耦合作用的結(jié)果。當間距比L*=1.5時，在某些折合流速(Ur=5.5)下，上游和下游圓柱的脫渦頻率不再相等，此時圓柱之間的耦合作用不再穩(wěn)定。此外，在某些折合流速下，三倍基頻占主導(dǎo)。