收稿日期：2022-06-25

基金項目：國家自然科學(xué)基金（520609010；51969009）

通信作者：羅竹梅（1979—），女，博士、教授、博士生導(dǎo)師，主要從事流固耦合及工程應(yīng)用理論方面的研究。112374102@qq.com

DOI：10.19912/j.0254-0096.tynxb.2022-0932 文章編號：0254-0096（2023）06-0045-08

摘 要：以剛性圓柱作為俘能結(jié)構(gòu)，通過強流固耦合數(shù)值方法模擬俘能結(jié)構(gòu)在低速水流環(huán)境下的三維渦激振動，對共振約化速度（[Ur=2～12]）范圍內(nèi)進(jìn)行俘能結(jié)構(gòu)的俘能大小及俘能效率分析。為獲得俘能結(jié)構(gòu)三維尾流對俘能大小的影響，采用熵產(chǎn)理論并結(jié)合尾流特性以捕捉俘能結(jié)構(gòu)尾渦的能量損失來源及分布，獲得不同振動分支下俘能大小與尾流耗散間的關(guān)系。結(jié)果表明，基于熵產(chǎn)理論的能量損失分析可準(zhǔn)確捕捉到俘能結(jié)構(gòu)的能量耗散趨勢。黏性熵產(chǎn)主要發(fā)生在圓柱表面，湍流熵發(fā)生在耗散率較大的尾流區(qū)。隨著約化速度的增加，上端分支內(nèi)，當(dāng)[Ur=6]時俘能結(jié)構(gòu)的俘能效率最大達(dá)到30.5%，俘能結(jié)構(gòu)俘獲的水動能增加，圓柱表面的黏性熵?fù)p耗和圓柱周圍及尾流區(qū)湍流熵產(chǎn)損失也增加，俘能結(jié)構(gòu)表面的黏性熵產(chǎn)損失和尾流區(qū)湍流熵產(chǎn)損失增大成為俘能效率下降的主要原因。

關(guān)鍵詞：海洋能；俘能結(jié)構(gòu)；渦激振動；流固耦合；能量耗散；熵產(chǎn)理論

中圖分類號：TK79；TK730.2 文獻(xiàn)標(biāo)志碼：A

0 引 言

實現(xiàn)碳達(dá)峰、碳中和目標(biāo)愿景，發(fā)展清潔可再生能源是重中之重。作為可再生能源之一的海洋能，自20世紀(jì)70年代開始，其開發(fā)就備受各沿海國家的重視。在海流流速穩(wěn)定的區(qū)域，可獲得穩(wěn)定可靠的海流能。目前人們對海流能的利用主要集中在潮流能的開發(fā)利用方面，即通過利用水下的渦輪發(fā)電設(shè)備可從有較大流速的潮流中獲取能量，但基于渦激振動從低速海流或河流中獲取能量方面的研究較少［1］。由于水下渦輪機受水流速度的限制，在低速水流中無法實現(xiàn)有效能量轉(zhuǎn)換，基于渦激振動驅(qū)動俘獲低速海流能或河流能成為清潔能源利用的一種新選擇［2］。

利用渦激振動從低速海流中俘獲能量的研究目前還停留在理論和實驗階段。相關(guān)文獻(xiàn)主要來源于Michael Bernitsas教授研究團(tuán)隊，該團(tuán)隊對從低速海流中俘獲海流能的能量轉(zhuǎn)換裝置VIVACE進(jìn)行了一系列實驗。Garcia［3］根據(jù)能量現(xiàn)象學(xué)思想建立出一個數(shù)學(xué)方程模型，將由Bernitsas團(tuán)隊對單個光滑的圓柱進(jìn)行的渦激振動的實驗數(shù)據(jù)代入方程，得到方程中的參數(shù)值。該數(shù)學(xué)模型還捕捉到渦激振動響應(yīng)的初始分支、上端分支和下端分支，但未捕捉到符合實驗結(jié)果的振幅，且模型中較小的能量變化都會對振幅產(chǎn)生較大影響。Bernitsas等［4］和Lee等［5］在高雷諾數(shù)和高阻尼條件下進(jìn)行渦激振動實驗，通過優(yōu)化俘能裝置的阻尼比，得到最大振幅比為1.78和最大能量轉(zhuǎn)換效率為30.8%以及平均能量轉(zhuǎn)換效率為22%。

利用數(shù)值方法進(jìn)行渦激振動俘能研究是Bernitsas團(tuán)隊的另一種研究手段，能輔助獲取俘能過程中流場變化信息與俘能參數(shù)和俘能大小間的關(guān)系，便于進(jìn)行試驗條件下難以實現(xiàn)的高雷諾數(shù)下的渦激振動響應(yīng)結(jié)果的研究。孫海等［6］采用二維非定常雷諾平均法（RANS），通過被動湍流控制對粗糙表面的單個彈性支撐圓柱的馳振現(xiàn)象進(jìn)行模擬，得到與實驗基本一致的響應(yīng)結(jié)果及渦脫落模式。丁琳等［7］采用OpenFOAM數(shù)值求解了兩個彈性支撐圓柱的橫向流激振動，采用被動湍流控制觀察到圓柱的馳振現(xiàn)象。振幅響應(yīng)除了初始分支、上端分支外，還有馳振分支。丁文俊等［8］采用相同的數(shù)值方法通過被動湍流控制方法研究了高雷諾數(shù)下的流激振動響應(yīng)，獲得較高的馳振振幅響應(yīng)。

僅有極少數(shù)的國外學(xué)者已對三維渦激振動進(jìn)行系統(tǒng)分析研究，Kondo［9］采用Petrov-Galerkin有限元法來求解三維的Navier-Stokes方程，模擬三維圓柱體的雙自由度渦激振動。Mittal［10］通過對質(zhì)量比[m*=10]和雷諾數(shù)[Re=1000]的圓柱體模型進(jìn)行三維數(shù)值模擬，重點研究約化速度[Ur]為4.2～6.2時上下分支間的滯后過渡。目前，對圓柱體中的渦激振動效應(yīng)進(jìn)行了大量系統(tǒng)的動力學(xué)實驗分析和計算數(shù)值研究，但大多數(shù)研究集中在抑制二維圓柱渦激振動響應(yīng)、減少工程破壞方面，對三維渦激振動俘能研究較少。

能量轉(zhuǎn)換過程中，不可避免地存在耗散效應(yīng)，局部熵產(chǎn)率（local entropy production rate，LEPR）正受到越來越多的關(guān)注。在整個流動區(qū)域上通過LEPR積分得到的產(chǎn)生熵與流動損失密切相關(guān)。Kock等［11］對管道流動進(jìn)行分析，指出熵產(chǎn)可作為一個有效的指標(biāo)，以確定傳熱的湍流過程的效率。此外，Herwig等［12］利用管道層流和紊流的實驗數(shù)據(jù)，驗證了熵產(chǎn)理論的正確性。通過對3種不同的二維粗糙度類型進(jìn)行模擬，計算任意形狀的壁面粗糙度所要求的等效沙粒粗糙度。在確定產(chǎn)生熵的基礎(chǔ)上，驗證了耗散模型適用于層流和湍流流動。

熵產(chǎn)理論普遍用于可壓縮流體中，主要用于評估熱力傳動的能量損失，而在不可壓縮過流部件中的應(yīng)用非常少。針對渦激振動現(xiàn)象，大多數(shù)研究者僅局限于二維渦激振動的水動力特性研究，真實低速海流環(huán)境下的三維渦激振動俘能研究較少，且基于熵產(chǎn)理論對渦激振動進(jìn)行能量捕捉及能量損失分析目前較少研究者涉及。為了研究渦激振動三維尾流效應(yīng)對俘能大小的影響、準(zhǔn)確捕捉圓柱四周及尾渦的能量變化規(guī)律，本研究采用流固耦合數(shù)值方法對均勻低速來流下的俘能結(jié)構(gòu)進(jìn)行三維渦激振動模擬，結(jié)合熵產(chǎn)理論，從能量損失的角度分析不同振動分支下俘能大小與尾流特性間的關(guān)系，尋找各分支下影響俘能能量耗散的主要原因。

1 渦激振動俘能模型

1.1 幾何模型

低速水流中渦激振動驅(qū)動發(fā)電的俘能結(jié)構(gòu)如圖1所示。整個俘能結(jié)構(gòu)為單圓柱構(gòu)成的質(zhì)量-彈簧-阻尼振動系統(tǒng)，結(jié)構(gòu)模型中振動圓柱直徑[D=0.1] m，長[L=0.5] m，質(zhì)量比[m*=2.4]，質(zhì)量阻尼比[m*ζ=0.013]。將圓柱置于均勻來流中，其在水中的固有頻率[fn=1.46 ]Hz，為模擬真實低速海流或河流環(huán)境，實驗流速[U=0.28～1.68] m/s，約化速度[Ur=2～12]。

1.2 數(shù)值計算模型

1.2.1 流體控制方程

均勻來流作用條件下俘能結(jié)構(gòu)的渦激振動屬流固耦合問題，選取RANS方程作為流體動力學(xué)控制方程，結(jié)合SST [k]-ω湍流模型對流場進(jìn)行基本理論方法求解，采用有限體積法對控制方程進(jìn)行離散，其流體控制方程為：

[?ui?xi=0] （1）

[?ρui?t+?ρuiuj?xj=-?p?xi+μ▽2ui-?ρui′uj′?xj] （2）

[-ρui′uj′=μt?ui?xj+?uj?xi-23ρktδij] （3）

式中：[ui]和[ui′]——對應(yīng)的[i]方向速度分量和速度變化；[xi]——[x]方向的笛卡爾坐標(biāo)系；[ρ]和[μ]——流體的密度和動力黏度；[t]——時間；[p]——壓力；[-ρui′uj′]——雷諾應(yīng)力；[μt]——湍流渦流黏度；[kt]——湍動能，從湍流振蕩過程產(chǎn)生的速度場中獲得；[δij]——Kronecker函數(shù)。

SST [k]-ω湍流模型當(dāng)中包括湍動能[k]和湍流耗散率[ω]的輸運方程。其中湍動能[k]的輸運方程為：

[?ρk?t+?ρujk?xj=P-β*ρωk+??xjμ+σkμi?k?xj] （4）

湍流耗散率ω的輸運方程為：

[?ρω?t+?ρujω?xj=γvtP-βρω2+" " " " " " ??xjμ+σωμt?ω?xj+21-F1ρσωω·?k?xj·?ω?xj] （5）

式中：[P]——雷諾應(yīng)力和平均運動變形率張量的二重標(biāo)量積；[β*]——熱膨脹系數(shù)；[σk]——湍動能Prandtl數(shù)；[γ、][β、][σω]——經(jīng)驗系數(shù)；[vt]——渦黏性系數(shù)；[F1]——混合函數(shù)，目的是使輸運方程從近壁面的[k]-ω模型過渡到外部的[k-ε]模型，其定義為：

[?=F1?1+1-F1?2] （6）

式中：[?1]——近壁面的[k]-ω模型，[γ1=β1/β*-][σω1κ2β*]；[?2]——外部的[k-ε]模型，[γ2=β2/β*][-σω2κ2β*]。其中，[σω1、][σω2、][β1、][β2]——在近壁面[k]-ω模型和外部[k-ε]模型中求解[k]和ω方程時的經(jīng)驗系數(shù)；[κ]——Von Karman常數(shù)。

1.2.2 固體控制方程

基于數(shù)值仿真方法，對自由表面下的圓柱振動方程進(jìn)行離散求解。固體運動方程采用四階龍格庫塔法，基于用戶自定義程序UDF求解每個時間步長的圓柱運動方程，其計算結(jié)果傳遞給流體方程，實現(xiàn)雙向流固耦合計算，該數(shù)值方法不僅計算效率較高，計算精度也較為滿意。四階龍格庫塔法求解圓柱運動方程如式（7）所示。

[fv，y=FfluidtM-2ξω0yt-ω02yt] （7）

采用四階龍格庫塔法將方程離散：

[vtn+1=vtn+Δt6×K1+2K2+2K3+K4] （8）

[ytn+1=ytn+Δt×vtn+Δt26K1+K2+K3] （9）

式中：[Δt]——時間步長；[K1]、[K2]、[K3]、[K4]——四階龍格庫塔的轉(zhuǎn)換因子。

如果確定了[tn]時刻的柱體結(jié)構(gòu)的位移[y（tn）]和速度[v（tn）]，可通過式（7）～式（9）來求解[tn+1]時刻的位移[y（tn+1）]。

1.3 計算域及網(wǎng)格劃分

為進(jìn)一步實現(xiàn)在真實來流條件作用下俘能結(jié)構(gòu)的三維渦激振動，建立出如圖2所示的流固耦合的計算模型。模型中俘能結(jié)構(gòu)為剛性單圓柱質(zhì)量-彈簧-阻尼系統(tǒng)，圓柱直徑為[D]，流體域長為[20D]，高為[5D]，寬為[10D]。圓柱中心距離上游進(jìn)口面為[5D]，距離出口面為[15D]。置于流體域中的俘能結(jié)構(gòu)在均勻來流激勵下，可捕捉低速水流環(huán)境下自由渦激振動的精細(xì)流場信息。

采用四邊形網(wǎng)格對流體的計算域進(jìn)行網(wǎng)格分布的精確劃分，沿水流流動方向，網(wǎng)格呈非均勻分布，為獲得圓柱表面及尾流精細(xì)的流場信息，對圓柱周圍流域的網(wǎng)格進(jìn)行加密，如圖3所示。流體域中的上下兩個邊界均為固壁界面，圓柱體與流體的交接界面處設(shè)置為流固耦合（fluid-structure-interaction，F(xiàn)SI）邊界。結(jié)構(gòu)表面采用無滑移邊界條件，即圓柱體表面的流體速度與圓柱振動速度相同，入口速度邊界設(shè)置為自由來流速度。下游邊界條件設(shè)定為壓力出口，參考壓力為零。

俘能結(jié)構(gòu)渦激振動產(chǎn)生較大位移，結(jié)構(gòu)的運動必然引起網(wǎng)格變形，本論文選取彈簧光順模型中的擴(kuò)散光順法（diffusion），通過求解式（10）的擴(kuò)散方程來設(shè)置網(wǎng)格節(jié)點：

[▽ι▽u=0] （10）

式中：[ι]——擴(kuò)散系數(shù)；[u]——網(wǎng)格的運動速度。

通過UDF自定義程序?qū)⑦\動圓柱設(shè)為剛性體，在流場中可以進(jìn)行整體運動，而不會發(fā)生彎曲扭轉(zhuǎn)和局部形變。

1.4 網(wǎng)格無關(guān)性驗證

為量化數(shù)值計算結(jié)果對網(wǎng)格密度的依賴性，對約化速度[Ur=5]時4種不同網(wǎng)格密度下俘能結(jié)構(gòu)的渦激振動進(jìn)行數(shù)值計算，結(jié)果見表1。由表1可看出，隨著圓柱表面徑向網(wǎng)格尺寸[Δr]的減少和總節(jié)點數(shù)[Nnode]的增加，振幅比[A*]和頻率比[f/fn]（[fn]為結(jié)構(gòu)在水中的固有頻率）結(jié)果逐漸趨于穩(wěn)定，采用第3種網(wǎng)格劃分方式既能滿足計算精度的要求，又能提高計算效率。

2 熵產(chǎn)理論

使用傳統(tǒng)的壓差評估方法無法準(zhǔn)確確定能量耗散的主要位置。圓柱的渦激振動過程中，由于工作介質(zhì)中的黏度和產(chǎn)生的雷諾應(yīng)力，機械能被不可逆地轉(zhuǎn)化為了內(nèi)能?；跓崃W(xué)第二定律，將熵產(chǎn)分析結(jié)合流體力學(xué)和傳熱理論來計算實際的能量損失，可準(zhǔn)確測量和捕捉能量損失。考慮到水的比熱容較高，三維圓柱渦激振動的數(shù)值計算中，認(rèn)為溫度為恒定，不考慮傳熱引起的熵產(chǎn)。在忽略熱傳導(dǎo)的情況下，邊界層中的黏性力將不可逆地轉(zhuǎn)化熱傳導(dǎo)流體的動能和勢能為內(nèi)能，從而引起熵的增加。

熵產(chǎn)主要包括黏度熵產(chǎn)、湍流熵產(chǎn)和溫度熵產(chǎn)。由于熵產(chǎn)是不可壓縮流體的狀態(tài)參數(shù)，引入熵產(chǎn)率來表示不可壓縮流體的輸運方程，如式（11）所示。

[ρ?s?t+u?s?x+v?s?y+w?s?z=-divqT+ΦT+ΦΘT2] （11）

式中：[s]——熵產(chǎn)；[u、v]和[w]——3個方向（[x、y]和[z]）的速度；[q]——熱流密度；[Φ]——黏性耗散產(chǎn)生的熵；[ΦΘ]——溫度變化產(chǎn)生的熵產(chǎn)。

不可壓縮流體的黏性耗散可通過式（12）計算。

[Φ=μ2?u?x2+?v?y2+?w?y2+" " " "?u?y+?v?x2+?u?z+?w?x2+?v?z+?w?y2] （12）

湍流中，局部熵產(chǎn)率（LEPR）由以下兩部分組成，其中一部分是時間平均運動所導(dǎo)致的熵產(chǎn)率即[S?D]，而另一部分是由速度波動導(dǎo)致的熵產(chǎn)率即[S?D′]。總熵產(chǎn)率計算公式為：

[S?D=S?D+S?D′] （13）

[S?D=2μeffT?u1?x12+?u2?x22+?u3?x32+μeffT" " " " " ?u2?x1+?u1?x22+?u3?x1+?u1?x32+?u2?x3+?u3?x22] （14）

[S?D′=2μeffT?u1′?x12+?u2′?x22+?u3′?x32+μeffT" " " " " ?u2′?x1+?u1′?x22+?u3′?x1+?u1′?x32+?u2′?x3+?u3′?x22] （15）

式中：[T]——溫度；[ui]——時間平均速度；[i]——笛卡爾坐標(biāo)系的3個方向；有效動態(tài)黏度用[μeff]表示，計算公式為：

[μeff=μ+μt] （16）

然而，在雷諾平均Navier-Stokes方法中，速度方向的[u1′，][u2′，][u3′]的分量不是直接計算得到的，而是通過求解[k]和ω方程間接得到的。根據(jù)Kock和Mathieu在[k-ω]湍流模式下的模擬計算，由速度波動引起的局部產(chǎn)生熵可近似計算如式（17）所示：

[S?D′=βρωkT] （17）

式中：[β]=0.09，湍流耗散速率[ε=βωk]，可對尾渦處的湍流耗散速率進(jìn)行分析。

因此，通過對整個流動區(qū)域的積分，可得到總產(chǎn)生熵。表達(dá)方式如式（18）～式（20）所示：

[SD=VS?DdV] （18）

[SD′=VS?D′dV] （19）

[SD=VS?DdV=SD+SD′] （20）

式中：[SD]——由時間平均運動引起的產(chǎn)生熵；[V]——流體區(qū)域的體積；[SD′]——由速度波動引起的產(chǎn)生熵；[SD]——全部產(chǎn)生熵。

3 結(jié)果分析

3.1 振動響應(yīng)結(jié)果

圖4為約化速度[Ur=2～12]范圍內(nèi)俘能結(jié)構(gòu)渦激振動響應(yīng)的振幅比[A*]和頻率比[f/fn]，采用與文獻(xiàn)［13］、文獻(xiàn)［14］數(shù)值實驗、文獻(xiàn)［15］物理實驗完全相同的低質(zhì)量比[m*=2.4]和低質(zhì)量阻尼比[m*ζ=0.013]，并將實驗?zāi)M結(jié)果與文獻(xiàn)結(jié)果加以比較。

從圖4可看出，本文采用的數(shù)值方法在各約化速度下的振幅比[A*]與文獻(xiàn)［15］的物理實驗結(jié)果基本吻合，但上端分支更寬，其中最大振幅比[A*]的值和物理實驗結(jié)果非常接近。當(dāng)頻率比[f/fn]在1.1～1.4之間時，鎖定現(xiàn)象發(fā)生，俘能結(jié)構(gòu)振幅較大，對應(yīng)于響應(yīng)的上端分支。此外，本數(shù)值實驗得到的鎖定范圍與文獻(xiàn)［14］、文獻(xiàn)［13］數(shù)值實驗結(jié)果基本一致，但振幅更高。鎖定范圍[Ur=4～10]內(nèi)，頻率比[f/fn]幾乎重疊。

3.2 三維尾流特性分析

渦激振動驅(qū)動的俘能結(jié)構(gòu)其尾流特性與其俘能結(jié)果有著必然聯(lián)系。捕捉俘能結(jié)構(gòu)的三維尾流特性可為結(jié)構(gòu)俘能中的能量分布及耗散規(guī)律提供流場分析依據(jù)。

圖5為不同約化速度下俘能結(jié)構(gòu)振動在振幅比最大位置時的三維尾渦脫落結(jié)構(gòu)，俘能裝置的尾渦結(jié)構(gòu)與來流速度和響應(yīng)幅值的相關(guān)性較大。[Ur=3]時，來流速度[U]僅為0.28 m/s，此時振幅響應(yīng)幅值較小，尾流區(qū)未出現(xiàn)明顯的脫落渦對，見圖5a。在渦激振動響應(yīng)鎖定范圍[Ur=4～10]內(nèi)，隨著約化速度的增加，流體與俘能結(jié)構(gòu)的耦合作用增強，圓柱結(jié)構(gòu)后出現(xiàn)明顯的脫落渦對。此時渦動力更大，響應(yīng)幅值較大，較大的響應(yīng)振幅導(dǎo)致出現(xiàn)更寬的尾流。由于響應(yīng)頻率[f]隨著約化速度增加而增加，尾流渦脫落頻率加快，尾流中脫落渦強度更強，相應(yīng)能量耗散也增加，如圖5b～圖5e所示。當(dāng)響應(yīng)進(jìn)入下端分支時，振動頻率偏離俘能結(jié)構(gòu)的固有頻率，響應(yīng)偏離鎖定范圍，導(dǎo)致渦激振動振幅減弱。在較大來流速度[U]和較低響應(yīng)振幅下，尾流區(qū)內(nèi)出現(xiàn)狹長的正負(fù)渦對，尾流寬度變窄，見圖5f。

截取不同約化速度下x-z平面內(nèi)的流場渦量圖，以清楚獲取尾渦模式隨約化速度的變化關(guān)系，如圖6所示。[Ur=3]時，響應(yīng)處于初始分支，此時出現(xiàn)一個經(jīng)典的2S尾流脫落模式，即在單個振動周期內(nèi)的尾流出現(xiàn)一對渦脫落，這與文獻(xiàn)［15］物理實驗中初始分支的渦脫落模式一致。上端分支內(nèi)，尾渦脫落模式從初始分支內(nèi)的2S模式變?yōu)?P模式，單個周期內(nèi)圓柱后出現(xiàn)一個較大渦對和一個較小渦對，隨著尾流向下游的發(fā)展，較小渦對逐漸耗散，見圖6b。上端分支內(nèi)隨著來流速度的增加，尾渦脫落頻率增強，流體結(jié)構(gòu)相互作用持續(xù)增強，在[Ur=9]時尾流中出現(xiàn)[P+T]尾渦模式，單個周期內(nèi)更多的渦對脫落。在下端分支中，渦位發(fā)生錯亂，展向渦排列非常不規(guī)則。尾渦模式變成P+S模式，圓柱上側(cè)脫離出兩個正渦，圓柱下側(cè)剛好脫落出一個完整的負(fù)渦。脫落的渦細(xì)長，占據(jù)下游較長的尾流區(qū)。

3.3 基于熵產(chǎn)理論的俘能分析

3.3.1 渦強度和耗散率分析

在圓柱渦激振動的過程中，產(chǎn)生的尾渦一部分轉(zhuǎn)換成圓柱渦激振動的升力，另一部分能量卻因為渦與渦之間的影響產(chǎn)生湍流能量耗散［16］。為了探究鎖定區(qū)內(nèi)俘能圓柱尾渦耗散規(guī)律，將[Ur=5]時俘能結(jié)構(gòu)尾流區(qū)捕捉到的4個脫落正渦中心[b1、b2、b3、b4]設(shè)為監(jiān)測點，提取4個監(jiān)測點的渦強度，計算由速度波動引起的局部產(chǎn)生熵[S?D′]中的[ε]，得到渦強度和尾渦耗散率[ε]隨時間變化的關(guān)系圖，如圖7所示。尾流渦主要由流向渦對控制，兩個相鄰渦對之間的距離約為1個圓柱直徑。從曲線圖中可清晰地看出，當(dāng)結(jié)構(gòu)產(chǎn)生較大振動幅值時，不斷有較強的渦對從圓柱上下表面脫落。由于振動位移隨時間發(fā)生周期性變化，第1對脫落的渦強度和耗散率隨時間也發(fā)生變化。當(dāng)脫落渦到監(jiān)測點[b1]處，此時渦強度和湍流耗散率分別為0.023和0.144，隨著脫落渦向下游不斷演變，渦強度不斷降低，湍流耗散率也越來越小，在[b4]點達(dá)到最低，此時渦強度和湍流耗散率分別降低了97%和93%。隨著時間的增加，越遠(yuǎn)離結(jié)構(gòu)的尾流區(qū)其尾渦強度和耗散率ε越小，渦耗散越快，流場受結(jié)構(gòu)的影響越小。

圖8顯示在4個約化速度（Ur=3、5、9、11）下第1個脫落正渦中心處的渦強度和耗散率三維曲線圖。結(jié)合圖4a和圖9可看出，渦強度及渦耗散率與不同約化速度下的振幅及俘能大小相關(guān)。當(dāng)俘能結(jié)構(gòu)處在較低振幅的初始分支（Ur=3）時，來流速度低，渦強度較小，為0.071，相應(yīng)的渦耗散率為0.063。隨著約化速度的增加，結(jié)構(gòu)受到的來流速度增加，在上端分支鎖定區(qū)間（Ur=5和Ur=9），流體結(jié)構(gòu)耦合作用持續(xù)增強，渦強度更大，俘獲的水動能[PVIV]增加。由于結(jié)構(gòu)振幅較大，尾流區(qū)流場受結(jié)構(gòu)干擾影響較大，相對于初始分支，尾渦耗散率增加了10倍。響應(yīng)進(jìn)入下端分支（Ur=11）時，結(jié)構(gòu)振動頻率偏離固有頻率，響應(yīng)進(jìn)入非鎖定區(qū)間，振幅降低，渦強度相應(yīng)也降低。此時大部分水動能都耗散在尾渦中，結(jié)構(gòu)從流體中俘獲的能量[PVIV]較少。

3.3.2 俘能結(jié)果分析

圓柱在均勻來流中俘獲的水動能采用式（21）［17］進(jìn)行計算。

[PVIV=1Tosc0ToscCtotaly2dt] （21）

式中：[Tosc]——單個周期振動時間；[Ctotal]——俘能結(jié)構(gòu)阻尼。

計算俘能結(jié)構(gòu)俘獲的的水動能[PVIV]和流體能量[Pfluid]的比值可獲得各約化速度下結(jié)構(gòu)的俘能效率[η]。對約化速度[Ur=2～12]范圍內(nèi)的俘能結(jié)構(gòu)渦激振動俘獲的水動能[PVIV]和俘能效率[η]進(jìn)行計算，得到如圖9所示的俘能結(jié)果和能量轉(zhuǎn)換效率隨約化速度變化關(guān)系圖，針對[Ur=3、][5、][9、][11]這4個典型約化速度（對應(yīng)圖中A、B、C、D），給出了圓周表面的渦黏性熵產(chǎn)損失圖。由圖9可知，初始分支中（[Ur=2～4]），來流速度較低，流體結(jié)構(gòu)干擾作用較弱，圓柱受到較小的升力作用，導(dǎo)致結(jié)構(gòu)從水流中獲取的水動能較少，為0.120 W。隨著約化速度增加，響應(yīng)進(jìn)入上端分支（[Ur=5～10]），此時鎖定區(qū)間振幅明顯增加，隨著渦動力的逐漸增強，流體結(jié)構(gòu)耦合作用增強，[Ur=10]時結(jié)構(gòu)從水流中獲取的能量達(dá)到最大為[6.486 W]。進(jìn)入下端分支后（[Ur=11～12]），此時來流速度更大，響應(yīng)處于非鎖定狀態(tài)，獲取的能量出現(xiàn)跳躍性減小，與上端分支相比降低了約80%，導(dǎo)致尾流場紊亂，圓柱表面因黏性耗散的熵產(chǎn)和尾流中因湍動能耗散的熵產(chǎn)都增加，水動能損失較大，俘能結(jié)構(gòu)僅獲取較少的水動能[PVIV]。對比能量轉(zhuǎn)換效率的曲線圖可知，隨著俘能結(jié)構(gòu)俘獲的水動能[PVIV]的增加，俘能效率[η]不斷增加。上端分支中，俘能最大效率達(dá)到30.5%。進(jìn)入下端分支，水流速度的增加使流體能量[Pfluid]隨之增加，由于[Pfluid]的增加速率高于俘能結(jié)構(gòu)在水流中俘獲水動能的增加速率，渦表面的黏性熵產(chǎn)損失也隨著流速的增加逐漸增大，導(dǎo)致能量轉(zhuǎn)換效率開始下降，圓柱渦激振動的黏性熵?fù)p耗增加，水動能損失增加，俘能結(jié)構(gòu)獲取的水動能減少，能量轉(zhuǎn)換效率也隨之減少。

3.3.3 熵產(chǎn)損失分布及能量耗散分析

能量耗散主要是因為渦動能和湍動能的耗散，這個過程常伴隨著速度梯度的巨大變化，從而導(dǎo)致尾流區(qū)不穩(wěn)定的流動。從圖8基于熵理論得出的尾渦強度和耗散率圖可知，水流撞擊圓柱表面時產(chǎn)生的能量損耗，黏性熵?fù)p耗和溫度熵?fù)p耗主要集中在圓柱與水的來流方向接觸面，尾流區(qū)有湍流熵?fù)p耗的產(chǎn)生。黏性熵產(chǎn)主要是由不可逆黏性耗散引起的，主要發(fā)生在圓柱靠近入口方向。由于均勻來流受圓柱附近渦的影響，改變了流動方向，形成渦流和沖擊損失，圓柱表面仍有小面積黏性耗散。與黏滯產(chǎn)生熵相比，湍流產(chǎn)生熵更大。

圖10將熵產(chǎn)理論用于評估三維渦激振動能量損失，給出了約化速度[Ur=3、][5、][9、][11]時流場的熵產(chǎn)分布圖，以分析各約化速度下渦激振動響應(yīng)與俘能結(jié)構(gòu)表面及尾流能量耗散間的關(guān)系。從圖10可看出，隨著約化速度的不斷增加，圓柱表面的黏性熵?fù)p耗和周圍的湍流熵?fù)p耗逐漸增加。當(dāng)俘能結(jié)構(gòu)響應(yīng)處于上端分支時（[Ur=5，Ur=9]），從圓柱后脫落的尾流旋轉(zhuǎn)角速度增大，旋渦強度隨切向速度的增大而增大。此時流動受到強烈干擾，黏度熵產(chǎn)和湍流熵產(chǎn)都比較大。能量轉(zhuǎn)換效率與圓柱渦激振動響應(yīng)振幅、頻率等息息相關(guān)，結(jié)合圖9能量轉(zhuǎn)換效率的曲線可知，[Ur=5]的時候能量轉(zhuǎn)換效率高，熵?fù)p耗小。由于尾渦處流體之間的動能交換，在俘能結(jié)構(gòu)的尾流區(qū)可觀察到數(shù)個湍流產(chǎn)生熵。

4 結(jié) 論

本文通過雙向流固耦合數(shù)值方法對均勻來流下三維圓柱渦激振動俘獲低速水流能進(jìn)行研究，利用熵產(chǎn)理論獲得了俘能結(jié)構(gòu)渦激振動響應(yīng)，俘能大小和俘能效率與能量損失分布及耗散間的關(guān)系，在本研究實驗參數(shù)范圍內(nèi)，得到如下主要結(jié)論：

1）三維效應(yīng)下，渦激振動的振幅頻率響應(yīng)結(jié)果與實驗結(jié)果接近，隨著約化速度的增加，俘能結(jié)構(gòu)尾流的流向渦對展向渦的控制作用更強。振幅越大，尾流區(qū)脫落渦橫向?qū)挾葘兏?，尾渦模式也從初始分支的2S模式變?yōu)樯隙朔种У?P和P+T模式，進(jìn)入下端分支，發(fā)展為P+S模式。

2）流體中的水動能其增加速率高于俘能結(jié)構(gòu)俘獲水動能的增加速率，隨著約化速度的增加，上端分支內(nèi)，當(dāng)[Ur=6]時俘能結(jié)構(gòu)的俘能效率最大達(dá)到30.5%，俘能結(jié)構(gòu)俘獲的水動能增加，而渦表面的黏性熵產(chǎn)損失增大，導(dǎo)致能量轉(zhuǎn)換效率逐漸下降。然而進(jìn)入下端分支后，圓柱渦激振動的湍流耗散熵增加，水動能損失增加較快，俘能結(jié)構(gòu)獲取的水動能驟減80%。

3）利用熵產(chǎn)理論準(zhǔn)確捕捉到能量損失的位置。隨著約化速度的增加，圓柱表面的黏性熵?fù)p耗和圓柱周圍及尾流區(qū)湍流熵?fù)p耗逐漸增加，且黏性熵產(chǎn)主要發(fā)生在圓柱表面，湍流熵發(fā)生在耗散率較大的尾流區(qū)。

引用熵產(chǎn)理論可有效預(yù)測和分析俘能結(jié)構(gòu)的能量損失，為俘能結(jié)構(gòu)設(shè)計和性能優(yōu)化提供參考，促進(jìn)俘能結(jié)構(gòu)在低速水流能中的應(yīng)用。

［參考文獻(xiàn)］

［1］ BARRERO-GIL" A，" ALONSO" G，" SANZ-ANDERS" A.Energy harvesting from transverse galloping［J］. Journal of sound and vibration， 2010， 329（14）： 2873-2883.

［2］ 羅竹梅， 聶聰， 郭濤. 渦激振動驅(qū)動的柱群結(jié)構(gòu)集中俘獲海流能研究［J］. 太陽能學(xué)報， 2021， 42（4）： 89-94.

LUO Z M， NIE C， GUO T. Study on centralized harvesting ocean current energy with column-group structure by VIV ［J］. Acta energiae solaris sinica， 2021， 42（4）： 89-94.

［3］ GARCIA E M H. Prediction by energy phenomenology for harnessing" "hydrokinetic" "energy" "using" "vortex-induced vibrations［D］. Michigan： University of Michigan， 2008.

［4］ BERNITSAS M M， BEN-SIMON Y， RAGHAVAN K， et al. The VIVACE converter：model tests at high damping and Reynolds number around 105［J］. Journal of offshore mechanics" "and" "arctic" "engineering，" "2009，" "131（1）： 011102.

［5］ LEE" J" H，" BERNITSAS" M" M." High-damping，" high-Reynolds VIV tests for energy harnessing using the VIVACE" converter［J］." Ocean" engineering，" 2011，" 38（16）： 1697-1712.

［6］ SUN" "H，" "BERNITSAS" "M" nbsp;M." "Bio-inspired" "adaptive damping in hydrokinetic energy harnessing using flow-induced oscillations［J］. Energy， 2019， 176： 940-960.

［7］ DING L， BERNITSAS M M， KIM E S. 2-D URANS vs. experiments of flow induced motions of two circular cylinders in tandem with passive turbulence control for 30， 000＜Re＜105，000［J］. Ocean engineering， 2013， 72： 429-440.

［8］ DING W J， SUN H， XU W H， et al. Numerical investigation on interactive FIO of two-tandem cylinders for hydrokinetic" "energy" "harnessing［J］. Ocean engineering， 2019， 187： 106215.

［9］ KONDO" "N." Three-dimensional" "computation" "for" "flow-induced vibrations in in-line and cross-flow directions of a circular" cylinder［J］. International" journal" for" numerical methods in fluids， 2012， 70（2）： 158-185.

［10］ MITTAL S. The critical mass phenomenon" in" vortex-induced vibration at low Re［J］. Journal of fluid mechanics， 2017， 820： 159-186.

［11］ KOCK F， HERWIG H. Local entropy production in turbulent shear flows： a high-Reynolds number model with wall functions［J］. International journal of heat and mass transfer， 2004， 47（10-11）： 2205-2215.

［12］ HERWIG H， GLOSS D， WENTERODT T. A new approach to understanding and modelling the influence of wall roughness on friction factors for pipe and channel flows［J］. Journal of fluid mechanics， 2008， 613： 35-53.

［13］ 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（4）： 449-466.

［14］ PAN Z Y， CUI W C， MIAO Q M. Numerical simulation of vortex-induced vibration of a circular cylinder at low mass-damping" "using" "RANS" "code［J］." Journal" of" "fluids" and structures， 2007， 23（1）： 23-37.

［15］ 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（7-8）： 813-851.

［16］ 譚俊哲， 王保振， 袁鵬， 等. 湍流強度對渦激振動潮流能轉(zhuǎn)換裝置振子響應(yīng)影響研究［J］. 太陽能學(xué)報， 2020， 41（10）： 20-26.

TAN J Z， WANG B Z， YUAN P， et al. Study on turbulence intensity influence on cylindrical oscillator response of VIV tidal energy conversion device［J］. Acta energiae solaris sinica， 2020， 41（10）： 20-26.

［17］ 羅竹梅， 張立翔. 影響從渦激振動中獲取能量的參數(shù)研究［J］. 振動與沖擊， 2014， 33（9）： 12-15.

LUO Z M， ZHANG L X. Influence of parameters on extracting energy from vortex-induced vibration［J］. Journal of vibration and shock， 2014， 33（9）： 12-15.

ENERGY CAPTURE ANALYSIS OF THREE-DIMENSIONAL VORTEX-INDUCED VIBRATION OF LOW-SPEED WATER FLOW BASED ON ENTROPY PRODUCTION THEORY

Li Jun1，Luo Zhumei1，Guo Tao2，Yang Tao1，Gao Suoming3

（1. School of Metallurgy and Energy Engineering， Kunming University of Science and Technology， Kunming 650093， China；

2. School of Engineering Mechanics， Kunming University of Science and Technology， Kunming 650500， China；

3. Kunming Branch of Shenzhen Water Planning and Design Co.， Ltd.， Kunming 650032， China）

Abstract：In this paper， the rigid cylinder is used as the energy-capture structure， and the three-dimensional vortex-induced vibration of the energy-captured system in the low-speed water flow environment is simulated by the strong fluid-solid coupling numerical method. The captured energy quantity and the energy capture efficiency of the capturing structure is analyzed" in the range of resonance reduction velocity（Ur=2-12）. To obtain the effect of the three-dimensional wake of the energy captured structure on the captured energy， the entropy production theory and the wake characteristics are used to capture the energy loss source and distribution of the wake vortex of the energy capture structure， and the relationship between capture energy and wake dissipation under different vibration branches are obtained. The results show that the energy loss analysis based on the entropy production theory can accurately capture the energy dissipation trend of the energy-harvesting structure. The viscous entropy production mainly occurs on the cylindrical surface， and the turbulent entropy occurs in the wake region， where the dissipation rate is more significant. With the increase of the reduced velocity， the kinetic energy of water captured by the energy capture structure in the upper branch increases， the viscous entropy loss on the cylinder surface and the turbulent entropy production loss around the cylinder and in the wake region also increases， and the viscous entropy production loss on the surface of the energy capture structure increases. The increase of turbulent entropy production loss in the wake region is the main reason for the decrease in energy capture efficiency.

Keywords：ocean energy； energy capture structure； vortex-induced vibration； fluid-structure interaction； energy dissipation； entropy generation theory