羅竹梅,張立翔
(1.昆明理工大學(xué) 能源與動力工程系,昆明 650093;2.昆明理工大學(xué) 工程力學(xué)系,昆明 650051)
目前對海流能的利用主要集中在潮流能利用方面,即利用水下渦輪發(fā)電裝置從具有較大流速的潮流中獲取能量,但基于渦激振動從海洋流或河流中獲取能量方面的研究很少[1-3]。非流線形結(jié)構(gòu)在水流作用下會在其兩側(cè)產(chǎn)生周期性的旋渦脫落,從而在結(jié)構(gòu)上形成與流向垂直的周期性外力,引發(fā)結(jié)構(gòu)橫向振動,即通常所說的渦激振動。渦激振動其共振發(fā)生在固有頻率附近一定范圍內(nèi),在較廣的流速和雷諾數(shù)范圍內(nèi),都能產(chǎn)生有效的振動,即使在低速的水流中也可獲得橫向的水動能。水下渦輪機在2.5~3.6 m/s甚至更高的流速下才可以工作,而全世界大量的海洋流其流速都低于1.5 m/s,一般河流其水流速度低于1.03 m/s。因此這一想法對于從低速河流和海流中獲取能量具有重要的現(xiàn)實意義。
利用渦激振動從低速水流中獲取能量的研究還停留在理論和實驗階段。Bernitsas[4-6]的團隊僅對單個彈性支撐剛性圓柱體的渦激振動進(jìn)行了研究,主要研究了自由表面、雷諾數(shù)、阻尼比及表面粗糙度的設(shè)計對圓柱體振幅和共振范圍的影響。圓柱體質(zhì)量比m*(m*為振動體質(zhì)量與排開流體質(zhì)量之比)、阻尼比ζ、質(zhì)量阻尼比m*ζ以及圓柱體在水中的固有頻率fn對從水流中獲取能量效率η的影響,目前相關(guān)文獻(xiàn)報道極少。受物理實驗條件的限制,本文通過流固充分耦合的數(shù)值方法,對二維彈性支撐圓柱體在均勻流中的橫向渦激振動進(jìn)行了模擬。并用該數(shù)值方法進(jìn)行了4種方案計算,以觀察m*、ζ、m*ζ和fn對η的影響。
流體控制方程用非定常不可壓縮雷諾平均N-S方程求解并用SST(Shear-Stress Transport)k-ω湍流模型封閉。為保證圓柱大位移運動時流體網(wǎng)格的質(zhì)量,采用任意拉格朗日歐拉(Arbitrary Lagrange-Euler,ALE)方法。ALE描述下流體的連續(xù)性方程和動量方程
(1)
(2)
(3)
采用有限元中的有限體積法對方程組進(jìn)行離散,流體區(qū)域采用非結(jié)構(gòu)四邊形網(wǎng)格,壓力-速度耦合方式采用PISO算法,空間導(dǎo)數(shù)采用二階精度的中心差分格式,壓力差值選擇標(biāo)準(zhǔn)方法,時間導(dǎo)數(shù)用二階組合時間積分法。
假定圓柱體被置于一個質(zhì)量-彈簧-阻尼系統(tǒng)中,如圖1所示。
圖1 圓柱體渦激振動模型
結(jié)構(gòu)控制方程采用有限元數(shù)值離散方法,離散后的非線性方程組寫成矩陣形式
(4)
流體和固體模型建好后,通過流固耦合求解器同時耦合求解,耦合邊界上應(yīng)用位移協(xié)調(diào)條件、力平衡條件以及速度平衡條件。每個時間步內(nèi)進(jìn)行多次迭代計算,直到滿足計算精度,流體和結(jié)構(gòu)才整體向前推進(jìn)計算。這一數(shù)值方法實現(xiàn)了流體和結(jié)構(gòu)充分耦合計算,可模擬真實的渦激振動現(xiàn)象。
為了驗證本數(shù)值方法的合理性,采用與Guilmineau等[7]數(shù)值實驗及Khalak等[8]的物理實驗完全相同的低質(zhì)量比m*=2.4和質(zhì)量阻尼比m*ζ=0.013。對振幅A、流速U和激勵頻率fex進(jìn)行無量綱化處理:A*=A/D,U*=U/fnD,f*=fex/fn,其中A*、U*和f*分別為圓柱體振幅比、約化速度和頻率比,實驗結(jié)果如圖2所示。
圖2 不同實驗中的振幅比A*
共振范圍內(nèi),圓柱體橫向位移響應(yīng)為定常態(tài)的諧波振動,呈正弦規(guī)律變化,且與升力之間存在一個相位差Φ,因此,可將單位長度圓柱體橫向位移y(t)和升力Fy表達(dá)成
y(t)=Asin(2πfext)
(5)
(6)
式中:Cy為升力系數(shù)幅值;D代表圓柱體直徑。
能量轉(zhuǎn)換效率η的大小可根據(jù)流體在一個振動周期內(nèi)對單位長度的圓柱體做功的大小Ps與單位長度內(nèi)流體的總能量Pf之比進(jìn)行計算,即:η=Ps/Pf。其中,單位時間內(nèi)單位長度流體的總能量
(7)
單位時間單位長度圓柱體從流體中獲取能量
(8)
結(jié)合式(5)~(8),采用如下無量綱形式表達(dá)能量獲取效率η
η=πA*CysinΦ(f*/U*)
(9)
約化速度U*的計算范圍為3~13,這一范圍包含了實驗中出現(xiàn)的整個共振范圍。用本文數(shù)值方法進(jìn)行了大量的數(shù)值計算,計算方案見表1。
表1 實驗方案
式(9)中各參數(shù)值可通過數(shù)值計算結(jié)果直接或間接得到,其中相位差Φ用過零鑒相法求得:在特定時刻,求解兩個信號通過零點的時間差Δt,根據(jù)參考信號的頻率ω,得出相位差Φ=180ωΔt/π,然后確定相位差Φ的符號。
為考察m*和ζ對橫向能量獲取效率η的影響,設(shè)計了兩組方案:第一組方案固定m*、fn,通過改變ζ以觀察ζ對η的影響;第二組方案則固定ζ、fn,通過改變m*以觀察m*對η的影響,見表1中方案1和方案2。
彈性支撐剛性圓柱體峰值振幅的大小取決于m*ζ,m*ζ越小,振幅峰值越大[8-10]。m*對渦激振動橫向位移的影響是非線性的,而結(jié)構(gòu)阻尼比ζ的影響則幾乎是線性的[11]。方案1的計算結(jié)果見圖3。
圖4中可以觀察到方案2下m*對η的影響趨勢:m*=2時,其能量獲取效率的峰值ηmax最小,為0.195,但從水流中有效獲取能量的約化速度最廣;m*=4時,其能量獲取效率的峰值ηmax最大,為0.262,相反,其從水流中有效獲取能量的約化速度范圍最窄。這說明m*不僅影響能量獲取效率的大小,同時也影響有效獲取能量的約化速度范圍。
圖3 方案1中A*和η隨U*的變化圖
從方案1的實驗結(jié)果不難發(fā)現(xiàn),當(dāng)m*ζ=0.21時,其ηmax出現(xiàn)最大值。為了觀察其它質(zhì)量阻尼比m*ζ下圓柱體從渦激振動中獲取能量的效率,增加了方案3,見表1。如圖5所示,當(dāng)m*=2,m*ζ=0.24時,其ηmax值遠(yuǎn)高于其它兩種阻尼比下的值,ηmax達(dá)到0.301,其它兩種情況的ηmax值相當(dāng),約為0.2。
圖4 方案2中η隨U*的變化圖
圖5 方案3中η隨U*的變化圖
圖6 m*ζ與η的關(guān)系圖
此外,在方案1的基礎(chǔ)上,增加了更小阻尼比ζ=0.01下的渦激振動響應(yīng)計算。將組合參數(shù)m*ζ值(0.03、0.1、0.15、0.2、0.21、0.24、0.3、0.4)及其對應(yīng)的ηmax值進(jìn)行整理,發(fā)現(xiàn)m*ζ控制ηmax值的大小,且存在一個最優(yōu)的m*ζ,使圓柱體從渦激振動中獲取能量的效率最大。對共振范圍內(nèi)各m*ζ值下的η求平均值,發(fā)現(xiàn)最大平均效率ηp也出現(xiàn)在m*ζ=0.24處,如圖6所示。因此,在某一雷諾數(shù)范圍內(nèi),選取合理的質(zhì)量阻尼比m*ζ,可使圓柱體從渦激振動中獲取能量的峰值效率ηmax及平均效率η最大。
圖7 方案4中A*隨U及η隨U*的變化圖
隨著fn的增加,其峰值振幅先增加后又減小,但其共振范圍卻向高流速和高雷諾數(shù)移動,這一點對從渦激振動中獲取能量很重要。與fn對峰值振幅的影響規(guī)律一樣,隨著fn的增加,其ηmax值也是先增加后減小,但有效獲取能量的約化速度范圍卻隨著fn的增加逐漸減小。從方案4的實驗結(jié)果圖7來看,適當(dāng)選取fn值,可使圓柱體從渦激振動中獲取能量的效率最大,且同時還可擁有較廣的速度范圍。
本文通過流固充分耦合數(shù)值方法模擬渦激振動,并用該數(shù)值方法對影響彈性支撐圓柱體從渦激振動中獲取能量的參數(shù)進(jìn)行了研究。設(shè)計了4種方案來分析m*、ζ、m*ζ和fn對獲取能量效率η的影響。通過對實驗數(shù)據(jù)進(jìn)行分析,得到如下結(jié)論:
(1)m*不僅影響η的大小,同時也影響有效獲取能量的約化速度范圍。當(dāng)m*不變時,存在一個最優(yōu)的ζ值使圓柱體從渦激振動中獲取能量的效率最大。因此,為了能在較大的速度范圍內(nèi)有效地從海洋流或河流中獲取能量,建議采用較低的m*,但為了得到更大的ηmax,需選擇一個合理的阻尼比ζ。
(2)m*ζ控制獲取能量的最大效率ηmax和平均效率ηp的大小,且存在一個最優(yōu)的m*ζ,可使ηmax和ηp最大。
(3)適當(dāng)選取fn,可使圓柱體從渦激振動中獲取能量的效率最大,同時還可使其應(yīng)用速度范圍較廣。
對從海洋流或河流中利用渦激振動獲取能量的振動系統(tǒng)設(shè)計,本實驗結(jié)果可提供參考。
[1]Barrero-Gil A,Alonso G,Sanz-Andres A.Energy harvesting from transverse galloping[J].Journal of Sound and Vibration,2010,329 (14):2873-2883.
[2]Bernitsas M M,Raghawan Y,Ben-Simon E M H,et al.VIVACE (vortex induced vibration for aquatic clean renewable energy from fluid flow[J].Journal of Offshore Mechanics and Arctic Engineering,2008,130: 041101.
[3]Yoshiki,Nishi.Power extraction from vortex-induced vibration of dual mass system[J].Journal of Sound and Vibration,2012:1-14.
[4]Bernitsas M M,Ben-Simon Y,Raghavan K,et al.The VIVACE converter: model tests at high damping and Reynolds number around105[J].Journal of Offshore Mechanics and Arctic Engineering,2009,131:011102.
[5]Raghavan K,Bernitsas M M.Experimental investigation of Reynolds number effect on vortex induced vibration of rigid circular cylinder on elastic supports[J].Ocean Engineering,2011,38:719-731.
[6]Lee J H,Bernitsas M M.High-damping,high-Reynolds VIV tests for energy harnessing using the VIVACE converter[J].Ocean Engineering 2011,38: 1697-1712.
[7]Guilmineau E,Queutey P.Numerical simulation of vortex-induced vibration of a circular cylinder with low mass-damping in a turbulent flow[J].Journal of Fluids and Structures,2004,19: 449-466.
[8]Khalak A,Williamson C H K.Motions,forces and mode transitions in vortex-induced vibration at low mass-damping[J].Journal of Fluids and Structures,1999,13: 813-851.
[9]范杰利,黃維平.細(xì)長立管兩向自由度渦激振動數(shù)值研究[J].振動與沖擊,2012,31(24):65-68.
FAN Jie-li,HUANG Wei-ping.Numerical simulation of 2-DOF vortex-induced vibration of a long riser[J].Journal of Vibration and Shock,2012,31(24):65-68.
[10]Govardhan R,Williamson C H K.Modes of vortex formation and frequency response of a freely vibrating cylinder[J].Journal of Fluids Mechanics,2000,420: 85-130.
[11]何長江,段忠東.二維圓柱渦激振動的數(shù)值模擬[J].海洋工程,2008,26:57-63.
HE Chang-jiang,DUAN Zhong-dong.Numerical simulation of two-dimensional cylindrical vortex-induced vibration[J].Journal of Ocean Engineering,2008,26:57-63.