半球形多向波聚合波道振蕩水柱氣室優(yōu)化設計?

杜小振， 常 恒， 趙 巖， 王 宇， 郝振富

(山東科技大學機械電子工程學院， 山東 青島 266590)

振蕩水柱(Oscillating water column 簡稱OWC)波浪能采集技術是利用起伏的海浪引起氣室內液面上下振蕩，再轉換為空氣動能驅動空氣透平系統(tǒng)發(fā)電輸出電能。主要研究方法為理論與數值水槽模擬和小型樣機實驗，結合入射波和海況條件設計最優(yōu)氣室結構。波浪作用于浮體的入射波分析：Vyzikas等[1]采用計算流體動力學(CFD)的物理場分析軟件(Open FOAM)的開源程序庫，研究了規(guī)則波和不規(guī)則波作用下數值波浪水槽內振蕩水柱運動。Kuo等[2]模擬分析了作用于OWC箱式防波堤上的不同入射波水動力特性，研究了波長比與空氣壓縮動力之間關系。Wang等[3]研究非線性波動和粘度對固定振蕩水柱裝置的水動力性能影響，確定流體非線性和粘性增加可提高空氣輸出動力。Medina-López E等[4]對較強的風暴氣候及海床變化進行數值模擬，結果表明OWC的效率隨海床特征尺寸增加而降低。氣室結構參數設計研究：寧德志等[5]基于域內源造波技術的時域高階邊界元方法， 建立了二維時域的非線性數值波浪水槽模型，研究了前墻入水深度、厚度與氣室寬度對結構共振和能量轉換效率的影響。Ahmed Elhanafi等[6-7]搭建1∶50比例的海上固定式OWC模型，采用試驗和仿真相結合方法研究了空氣壓縮與氣室波浪俘獲能力關系，同時分析了最大效率時入射波高與前墻吃水深度的關系。胡杭輝等[8]基于勢流理論建立了二維完全非線性波作用下OWC裝置的數值模型，總結了波高、裝置吃水深度、氣室寬度和墻體厚度對能量轉換效率的影響。He等[9]將雙氣室OWC結構嵌入浮箱式防波堤擴大能量捕獲帶寬，研究波浪周期、吃水深度和腔室布置等對輸出功率影響。Ahmed Elhanafi等[10]研究前墻吃水深度、氣室寬度和動力輸出(POWER-TAKE-OFF，PTO)阻尼對雙氣室OWC裝置波浪能提取效率影響。海洋浮標多為圓柱形或圓錐形結構，波浪方向隨風向、洋流、潮汐等不斷變化，浮標式OWC波浪能捕獲結構需具備多向迎波，為此，周宇等[11]提出可以吸收各個方向海浪的“蘑菇型”氣室結構，討論了氣室寬度、外墻厚度和吃水深度對波面幅值變化和波能俘獲系數影響。寧德志等[12]研究圓柱形氣室內的自由表面振蕩模式，較低吃水深度和增加腔室寬度可以增加共振頻率擴大振蕩頻率帶寬。目前，振蕩水柱波浪能采集氣室研究了入射條件、尺寸結構優(yōu)化計算和仿真分析，結構以固定單一方向波浪能采集結構設計為主。為了提高浮標振蕩水柱是波浪能捕獲效率，本文提出半球形多向聚合波道迎波振蕩水柱氣室結構模型，研究提升遠海點式吸收波浪發(fā)電效率方法。

1 波浪運動數值模擬計算控制方程

振蕩水柱氣室處于三維數值波浪水槽內，假定流體為黏性、不可壓縮、自由流動，則數值波浪運動基本控制方程由連續(xù)性方程和動量方程組成[13]。

連續(xù)性方程：

(1)

動量方程：

(2)

(3)

(4)

式中： 流體密度u、v和w分別代表沿x、y和z方向的速度分量；μ為流體動力粘性系數；p為流體內部壓力；g為重力加速度；t為時間。

采用流體體積法(Volume of Fluid，VOF)[14]追蹤數值波浪水槽內氣液分層兩相流的兩個流體交界面，以確定每個控制容積單元中流體(空氣和水)的體積分數。體積分數連續(xù)方程：

(5)

式中：U是速度矢量；α是水的體積分數，當α滿足0<α<1時，控制容積單元中存在流體交界面。基于體積分數計算水氣混合密度:

ρ=αρw+(1-α)ρα。

(6)

式中：ρw是海水密度；ρα是空氣密度。

2 三維數值水槽仿真分析

2.1 半球多向迎波振蕩水柱氣室結構設計

振蕩水柱氣室優(yōu)化設計可有效提高波浪能轉化效率。以傳統(tǒng)的振蕩水柱波浪能采集氣室設計原理為基礎，提出半球形聚合波道多向迎波離岸式振蕩水柱波能轉換氣室，氣室上端形狀能有效解決強海風作用下的發(fā)電裝置的穩(wěn)定性，氣室內部增設后壁如圖1(a)，三種不同開口角度(β=120°、90°和72°)波道后壁模型如圖1(b)，探索多角度后壁的多向波采集技術，增加氣室內壓強、出氣口速度等提高發(fā)電效率。

2.2 數值水槽和氣室三維建模

利用NX.UG三維建模軟件建立1∶1比例的振蕩水柱氣室及三維水槽模型，將氣室模型置于數值水槽內部，計算區(qū)域模型如圖2。

圖1 半球形聚合波道多向迎波振蕩水柱氣室結構模型

圖2 三維數值水槽和氣室模型

水槽長L0=42 m，寬W0=6 m，高H0=8 m(分為空氣相4 m和水相4 m)。充分考慮造波區(qū)、消波區(qū)長度及水槽側壁效應，將水槽寬度設置為OWC氣室直徑的2倍，氣室中心距造波邊界21 m，距水槽兩側壁面各3 m，氣室底端距水槽底部0.5 m。在滿足數值造波功能及計算準確性的基礎上，合理設計模型的幾何尺寸可以減少計算量，實現快速仿真分析，氣室和數值水槽具體參數如表1所示。

三維數值水槽和OWC氣室模型ICEM CFD結構化網格劃分過程采用混合網格法，將水槽分為左側造波區(qū)，中部氣室區(qū)，右側消波區(qū)；其中邊界條件設置如圖3(a)所示，水槽左側造波邊界設置為速度入口，右側消波邊界設置為壓力出口；兩側面和底部邊界均設置為固壁表面；頂部邊界受恒定大氣壓力的影響，設置為壓力入口；左側造波區(qū)和右側消波區(qū)為形狀規(guī)則的長方體，故對其進行結構化網格劃分，并將波浪運動區(qū)域做網格加密處理，如圖3(b)所示；中部氣室區(qū)需要裝配氣室，采用非結構化網格劃分生成高質量的網格,如圖3(c)所示。

表1 三維數值水槽及氣室結構參數

圖3 三維數值水槽和氣室網格劃分

2.3 網格獨立性驗證

網格獨立性驗證，恰當劃分網格，既保證數值模擬仿真精度又降低計算時間。模擬規(guī)則波正向入射，在Fluent中監(jiān)測振蕩水柱氣室頂部壓強進行收斂測試以確定合適的ICEM網格。因三維數值水槽的長度遠大于其高度和寬度，文中僅改變x軸方向結構化網格數量，驗證沿波長方向網格大小對仿真結果的影響，介紹四種接近收斂的結構化網格測試方案如表2所示。從方案1到方案4，氣室造波區(qū)和消波區(qū)的結構化網格沿x軸方向逐漸加密，氣室工作區(qū)的非結構網格設置全局網格尺寸系數(Scale factor)和最大網格尺寸(Max element)分別為2.0和0.1。

表2 結構化網格數量

四種不同網格密度計算所得的振蕩水柱氣室頂部壓強曲線如圖4，不同網格密度的壓強分析結果保持一致；網格相對密集的方案2～4計算所得的壓強峰值基本相同；網格數量較少的方案1計算所得壓強略小于其它方案，平均誤差約為8.9%且僅出現在波峰和波谷部分。因此，方案1網格較為稀疏，雖然能夠模擬趨勢，但計算精度較低，比較方案2與方案3和4，則保證精度條件下，本文數值水槽及振蕩水柱氣室仿真模擬研究均采用方案2劃分網格。

圖4 網格獨立性驗證

2.4 數值水槽波形驗證

基于VOF流體體積法及N-S控制方程建立三維仿真分析模型，波浪模擬采用流線彎曲、漩渦模擬精確的RNG模型，選擇明渠造波法和多孔介質消波法分別實現造波和消波。根據渤海、萊州灣口等近海計算工況，仿真分析過程定義計算域內規(guī)則波的平均水深4 m、波高0.8 m、入射波周期3.5 s、波長17.179 5 m、湍流強度0.12%。利用分區(qū)標記功能(Region功能)劃分空氣和水的計算域，定義空氣相為第一相，水相為第二相；利用Pach功能定義水的體積分數1。分別在距水槽左側造波邊界1、6、12和18 m處設置監(jiān)測面，記錄水槽和氣室內波浪高度變化，位置如圖5，設置仿真時間步長0.01 s，迭代步數4 000步。

圖5 監(jiān)測點和監(jiān)測面示意圖

在三維數值水槽中模擬二階斯托克斯波，可通過設置在數值模型中的監(jiān)測面1～4對水槽不同位置的波浪上升時間、波峰和波谷形狀等波浪特性進行監(jiān)控，對波浪穩(wěn)定性、波高沿程衰減量及周期誤差進行分析。對于波高H=0.8 m、波周期T=3.5 s、波長=17.179 5 m的非線性規(guī)則波，數值水槽4個位置的自由表面波高分別如圖6(a)～(d)所示，將仿真模擬獲得的數值波與二階Stokes波的理論值進行對比，驗證放置振蕩水柱氣室的數值水槽波形準確性。

對比同一監(jiān)測面處數值波與二階Stokes波的波形曲線，可知數值模擬波的自由表面高度隨時間發(fā)展逐漸升高并接近于理論值，二者波浪周期基本吻合且波形均呈余弦變化。對比不同監(jiān)測面處的數值波波形，可知監(jiān)測面距離左側邊界越遠，波形穩(wěn)定所需時間越長，波浪高度沿程衰減越大。分別取16～30 s內x=1 m和x=18 m處波高的平均值，x=18 m處波浪高度較x=1 m處波浪幅值減小約14.33%，由于仿真模擬過程計入流體粘性作用和非線性影響，使得波浪在傳播過程中存在數值耗散的衰減情況。綜上，數值水槽中波浪高度的仿真結果和理論結果吻合，三維數值水槽模擬波高準確，可用于振蕩水柱氣室的仿真分析。

圖6 數值波和理論波對比曲線

2.5 波浪與氣室相互作用分析

數值模擬有限水深條件下振蕩水柱氣室與正向入射規(guī)則波浪的相互作用，將計算域內的空氣相與水相分區(qū)并定義水的體積分數為1，得到數值水槽的體積分數云圖如圖7。上部藍色區(qū)域為空氣相，下部紅色區(qū)域為水相，水氣交界面清晰的展示了自由振動波面。觀察自由表面的波形變化，左側造波區(qū)形成周期性規(guī)則變化的波浪形狀，右側消波區(qū)波浪趨于平穩(wěn)，數值水槽實現了有效的數值造波與消波。

圖7 數值水槽體積分數云圖

為分析波浪與氣室的相互作用，建立水槽中截面觀察氣室內外的自由液面變化，隨時間變化每1/7波周期的體積分數云圖如圖8所示。觀察圖8(a)～(c)，隨入射波波峰運動到中部氣室區(qū)，氣室外波面逐漸升高并向氣室內傳遞，氣室內波面隨之升高，形成上升水柱；觀察圖8(d)～(f)，隨入射波波峰向波谷過度，氣室外波面逐漸下降，氣室內波面也隨之下降，形成下降水柱。但因氣室外壁對波浪傳播具有阻礙作用，氣室內波面上升(或下降)的時間均晚于氣室外波面上升(或下降)的時間。由上可知，氣室內部水氣交界面隨波浪周期變化產生振蕩水柱，實現了仿真過程中波浪與氣室的相互作用。

當t=75 s時，數值水槽及振蕩水柱氣室內的流場渦量變化如圖9(a)，設定參考渦量數值為0 s-1，整個計算域內流場的渦量均為正數，其中波峰、波谷、氣室外壁及氣室內部的渦量較大。因多相流仿真包含多種流體介質，所以分開討論水域和空氣域中的流場渦量變化，如圖9(b)、(c)所示。水域中氣室外壁底端渦量較大，這是因為氣室外壁相對于水流的剪切運動產生大量渦量，在流體粘性作用下逐漸耗散；空氣域中氣室出氣口處產生渦量較多，最大值為76.93 s-1，振蕩水柱擠壓空氣發(fā)生旋轉運動形成渦流，一部分通過出氣口向氣室外運動和耗散，另一部分在氣室內部聚集成渦量團，渦量大小清楚顯示出了在波浪與氣室相互作用下流體速度場的旋度特征。

圖8 氣室每1/7波周期的體積分數云圖

圖9 數值水槽和振蕩水柱氣室渦量云圖

3 氣室模擬結果與討論

在規(guī)則波條件下，數值波浪水槽模擬波浪作用于半球形振蕩水柱氣室，形成空氣壓強并且出氣口速度隨入射波進入氣室逐漸增大，對表1所示的三維數值水槽及氣室仿真參數計算結果顯示，初次波浪進入氣室5.0 s后氣壓趨于穩(wěn)定的周期性變載荷狀態(tài)。因此氣室特性分析以此時間段后產生的壓強和速度為主。研究氣室后壁對半球形氣室波浪俘獲效率的影響，優(yōu)化波道開口角度實現離岸式振蕩水柱氣室多向聚波采集技術。

3.1 氣室氣流運動特性分析

圖7按時間順序給出氣室出氣口的速度矢量圖，顯示了一個周期內氣室出氣口的速度大小和氣流走向：圖10(a)～(c)，氣室內氣體因入射波波峰進入氣室受擠壓排出，排氣速度先增大后減小，隨入射波波峰向波谷過度，氣室由排氣狀態(tài)轉換為吸氣狀態(tài)，氣室內氣流紊亂；圖10(d)～(f)，氣室內自由液面隨入射波波谷進入下降，氣室因內、外壓力差變化開始吸氣，吸氣速度先增大后減小，最后因入射波波谷向波峰過度，氣室內氣流再次紊亂。

3.2 增設聚波后壁對氣室性能影響

半球形氣室中心迎波后壁對氣室壓強和出氣口速度變化影響曲線如圖11。在波浪穩(wěn)定后不同氣室結構的壓強變化規(guī)律基本相同，壓強幅值周期與入射波周期吻合，氣室增設聚波后壁產生的空氣壓強顯著高于無后壁氣室結構，增設90°后壁時，波浪進入氣室產生最大正壓約為1 550 Pa，產生最大負壓約為1 900 Pa。氣室聚波后壁的加入有利于提高半球形氣室出氣口的速度，增設氣室后壁后最大出氣口速度可以達到53 m/s，這是因為增設氣室后壁使波道口捕獲的波浪增加，提高氣室聚波能力，氣室內波面振蕩加強提高了空氣壓強和出氣口速度。

圖10 單周期內不同時間點氣室出氣口速度矢量圖

圖11 后壁對氣室壓強和出氣口速度的影響

3.3 波道開口角度對氣室性能的影響

假定各氣室模型波道開口均朝向來波方向，氣室數值變化規(guī)律基本相同；三種不同波道開口角度(β=120°、β=90°、β=72°)氣室模型的壓強對比曲線如圖12，當β=120°時氣室壓強明顯高于其它氣室，波浪進入氣室產生最大正壓約為1 850 Pa，達到波谷時則產生最大負壓約為2 300 Pa。氣室出氣口速度對比曲線如圖13，氣室正壓峰值所對應的氣室出氣口最大速度為V120°>V90°>V72°，氣室負壓峰值所對應的氣室出氣口最大速度基本相同，減小氣室波道開口角度將引起氣室正壓峰值對應的出氣口速度降低，相鄰速度峰值落差較大，不利于能量的收集。隨波道開口角度減小，入射波集聚反射增加，波能捕獲能力減弱。對比分析，氣室聚波后壁開口β=120°時氣室出氣口氣壓速度穩(wěn)定，有利于波浪能捕獲和轉換。

圖12 氣室壓強與l迎波角度關系

圖13 迎波角度氣室出氣口速度對比圖

3.4 前壁形狀對氣室性能的影響

由振蕩水柱氣室渦量云圖可知，前壁相對于水流的剪切運動產生大量渦量，部分能量在流體粘性作用下逐漸耗散。前壁作為波浪能轉化裝置的能量入口，合理的前壁厚度和入水深度可以有效減小能量損失，提高氣室的波浪采集效率。選擇最優(yōu)后壁開口角度的基礎上，選取4種前壁入水深度和3種前壁厚度兩兩組合，建立12種不同形狀的前壁，如表3所示。在三維數值波浪水槽中對12種氣室結構進行模擬仿真，對比研究不同前壁入水深度和不同前壁厚度對半球形振蕩水柱氣室工作性能的影響。

在相同規(guī)則波條件下分別對12種不同氣室結構進行模擬仿真，得到三種不同前壁厚度下的氣室壓強和出氣口速度變化曲線分別如圖14～16所示。當振蕩水柱氣室壁厚分別為0.1、0.3 和0.5 m時，不同前壁入水深度下的氣室壓強和出氣口速度曲線變化規(guī)律基本相同；氣室壁厚保持一定，壓強和速度均隨入水深度的增加而不斷減小，即三種壁厚的氣室均在入水深度h0=0.5 m時獲得最大頂部壓強和出氣口速度。例如，當振蕩水柱氣室的前壁厚度b=0.5 m保持不變時，4種不同入水深度下的氣室頂部壓強依次為Ph0=0.5 m>Ph0=1 m>Ph0=1.5 m>Ph0=2 m，出氣口速度依次為Vh0=0.5 m>Vh0=1 m>Vh0=1.5 m>Vh0=2 m，入水深度h0=0.5 m的振蕩水柱氣室具有最大頂部壓強1 804.50 Pa及最大出氣口速度59.03 m/s。由上可知，在一定范圍內減小前壁入水深度可有效減小波浪觸壁反射的現象，提高氣室頂部壓強和出氣口速度，增強振蕩水柱氣室的波浪俘獲能力。

表3 氣室前壁尺寸參數

當t=10.5～17.5 s時，3種前壁厚度及4種入水深度下的氣室壓強和出氣口速度的最大值如圖17所示。當前壁入水深度分別設定為0.5、1、1.5、2 m 4種情況時，相同波浪條件下的氣室壓強和速度均隨氣室壁厚的增加而增加，在壁厚為0.5 m時獲得最大值。例如，當入水深度h0=0.5 m時，3種壁厚(b=0.1、0.3、0.5 m)氣室獲得最大頂部壓強分別為1 697.88、1 798.21、1 804.50 Pa，最大出氣口速度分別為50.42、57.10、59.03 m/s。由上可知，當前壁入水深度一定時，隨氣室壁厚的逐漸增加，仿真所得頂部壓強和出氣口速度均不斷增大。

圖14 前壁入水深度對出氣口空氣壓強和速度的影響(b=0.1 m)

圖15 前壁入水深度對壓強和速度的影響(b=0.3 m)

圖16 前壁入水深度對壓強和速度的影響(b=0.5 m)

由以上分析可知，增大氣室前壁厚度或減小前壁入水深度，可有效減小波浪觸壁反射形成的能量耗散，有利于提高振蕩水柱氣室對波浪能的捕獲能力。對比12種前壁結構，當壁厚b=0.5 m、前壁入水深度h0=0.5 m時半球形氣室結構最優(yōu)，仿真獲得頂部壓強和出氣口速度的最大值分別為1 804.50 Pa和59.03 m/s。

4 結論

在規(guī)則波正向入射條件下對半球形振蕩水柱氣室進行數值模擬，未考慮復雜的海洋環(huán)境?；赩OF流體體積法及N-S控制方程建立半球形振蕩水柱波浪能采集三維數值水槽的模擬分析模型，利用數值模擬方法研究了半球形多向迎波振蕩水柱氣室能量俘獲特點，依據分析結果分析了增設氣室后壁、改變波道開口角度、優(yōu)化前壁形狀對半球形氣室性能的影響。

(1)增設氣室后壁的半球形氣室的壓強和出氣口速度均有提高，氣室對波浪能的俘獲能力得到提升。

(2)減小波道開口角度將降低氣室正壓峰值和出氣口速度，波道開口β=120°時氣室壓強最大、出氣口速度穩(wěn)定。

(3)增大氣室前壁厚度或減小前壁入水深度，波浪觸壁反射形成的能量耗散減小，振蕩水柱氣室對波浪能的捕獲能力得到提升。