養(yǎng)殖工船養(yǎng)艙流場數(shù)值仿真與特性分析

王銀濤，李志松,2，郭曉宇,2,，崔銘超

(1上海交通大學(xué)船舶海洋與建筑工程學(xué)院，上海 200240；2上海交通大學(xué)水動(dòng)力學(xué)教育部重點(diǎn)實(shí)驗(yàn)室，上海 200240；3 中國水產(chǎn)科學(xué)研究院漁業(yè)機(jī)械儀器研究所，上海 200092)

隨著社會(huì)的發(fā)展，以豬、牛為代表的“紅肉”消費(fèi)逐漸將被以雞、魚為代表的“白肉”代替，聯(lián)合國糧農(nóng)組織(FAO)的統(tǒng)計(jì)數(shù)據(jù)顯示[1]：中國的水產(chǎn)品產(chǎn)量逐年遞增，在2020年已經(jīng)達(dá)到了8 392.9萬t。但是由于近岸和陸基養(yǎng)殖的空間限制以及環(huán)境污染，漁業(yè)養(yǎng)殖正在歷經(jīng)從內(nèi)陸、沿海走向深遠(yuǎn)海養(yǎng)殖的轉(zhuǎn)變過程[2]。專業(yè)化的大型養(yǎng)殖工船集養(yǎng)殖、加工、儲存、物流于一體，可有效利用深遠(yuǎn)海適宜海域進(jìn)行長期游弋式養(yǎng)殖，并能躲避臺風(fēng)、赤潮等自然災(zāi)害，被稱為“移動(dòng)的海洋牧場”，工船養(yǎng)殖是發(fā)展深遠(yuǎn)海養(yǎng)殖的一個(gè)重要的方向[3-4]。

從養(yǎng)殖工船的概念提出到全球首艘10萬噸級深遠(yuǎn)海養(yǎng)殖工船的建設(shè)完成并投入運(yùn)行，國內(nèi)研究人員主要在總體規(guī)劃和設(shè)計(jì)等方面開展了研究[5-9]。養(yǎng)殖工船內(nèi)的養(yǎng)殖水艙在波浪激勵(lì)下產(chǎn)生晃蕩，Guo等[10]、Cui等[11]研究了極端工作環(huán)境下養(yǎng)艙內(nèi)流場特性，由于進(jìn)水射流在劇烈晃蕩下并不明顯，因此忽略了進(jìn)出水的影響。養(yǎng)殖工船多數(shù)時(shí)間在海況溫和的環(huán)境下作業(yè)，進(jìn)出水形成的推流作用明顯，研究“國信1號”實(shí)體養(yǎng)艙在進(jìn)出水影響下流場特性具有現(xiàn)實(shí)意義。

目前，養(yǎng)艙進(jìn)出水流場特性研究主要集中在陸基養(yǎng)殖，研究指出養(yǎng)艙形狀與進(jìn)出水條件都會(huì)對流場特性產(chǎn)生影響[12-16]。與工船養(yǎng)艙不同，陸基養(yǎng)艙內(nèi)豎有直立柱提供進(jìn)水，且通常水深較淺。且由于進(jìn)出水推流與自由面影響，養(yǎng)艙內(nèi)會(huì)產(chǎn)生漩渦流，進(jìn)水形式對漩渦流場影響以及流場適漁性目前未有相關(guān)研究。

漩渦流機(jī)理復(fù)雜，研究者基于對旋渦的認(rèn)識提出了理論模型[17-19]，但是理論分析僅適用于簡單的幾何情況。也有研究者基于不同的應(yīng)用背景開展了模型試驗(yàn)，研究結(jié)果表明漩渦形成的臨界深度和弗勞德數(shù)、雷諾數(shù)、韋伯?dāng)?shù)存在相關(guān)性[20-23]。隨著計(jì)算機(jī)技術(shù)發(fā)展，數(shù)值模擬成為漩渦流場研究的有效手段，通過發(fā)展數(shù)值模型研究自由表面渦旋流場的分布[24-26]，從數(shù)值角度證實(shí)了自由表面漩渦是類蘭金渦。漩渦流的模擬依賴湍流模型，采用大渦模擬方法可得到不同紊流強(qiáng)度以及邊界條件下進(jìn)水非穩(wěn)態(tài)漩渦運(yùn)動(dòng)[27]。由于現(xiàn)有的湍流模型耗散過強(qiáng)，一些湍流模型不能有效地捕捉自由表面渦[28]，目前研究表明大渦模擬(LES)[29]和SST(Shear Stress Transport)模型[30-31]可用于自由表面渦的模擬。開源軟件OpenFOAM的InterFoam求解器被廣泛應(yīng)用于氣液兩相流的求解，一些研究使用Openfoam模擬自由表面渦[32]，通過對比試驗(yàn)和理論，驗(yàn)證了其在漩渦求解方面的適用性。

以上關(guān)于自由表面漩渦流的研究主要還是基于較淺的水(水深和出水口直徑比較小)，研究的目的主要在于如何防止漩渦的產(chǎn)生，關(guān)于漩渦流的研究暫未涉及養(yǎng)殖工船養(yǎng)艙適漁性。為優(yōu)化養(yǎng)艙結(jié)構(gòu)與進(jìn)出水系統(tǒng)，養(yǎng)艙漩渦流場亟待研究。本研究基于OpenFOAM采用DES湍流模型對大型養(yǎng)殖工船養(yǎng)艙中的漩渦流進(jìn)行研究，分析流場速度特性以及適漁性，并通過改變?nèi)胨诘倪吔鐥l件來改變進(jìn)水流量和角度，分析流量和角度對漩渦流場的影響，為工程應(yīng)用提供參考。

1 數(shù)學(xué)模型

1.1 控制方程

不可壓縮流的連續(xù)性方程和動(dòng)量方程分別為：

(1)

(2)

式中：u為速度，m/s；ρ為密度，kg/m3；μ為運(yùn)動(dòng)粘性系數(shù)，kg/(m·s);F為體力，kg/(m2·s2)；p為壓強(qiáng)，kg/(m·s2)。

湍流模型采用kOmegaSSTDES模型，對于 kOmegaSSTDES模型，湍動(dòng)能k和比湍流耗散率ω的輸運(yùn)方程分別表示為[33-35]：

(3)

(4)

其中：

(5)

(6)

(7)

式中：F1為SST模型的邊界層識別函數(shù)；Γk、Γω分別為湍動(dòng)能和比湍流耗散率的有效擴(kuò)散系數(shù)；Gk、Gω分別為湍動(dòng)能和比湍流耗散率的湍流生成項(xiàng)；Ek、Eω分別為湍動(dòng)能和比湍流耗散率的耗散項(xiàng)。各項(xiàng)表達(dá)式和參數(shù)參照文獻(xiàn)[35]。

方程(3)中湍動(dòng)能k的耗散項(xiàng)Ek=-ρβ*ωkFDES，其中FDES為SSTDES的開關(guān)函數(shù)，可以表示為：

FDES=max{lk-ω/(CDESΔ),1}

(8)

lk-ω=k1/2/(β*ω)

(9)

CDES=0.61,Δ=max(δx,δy,δz)

(10)

式中：lk-ω為k-ωSST模型中長度尺度；CDES為模型參數(shù)。

1.2 多相流模型

OpenFOAM中interFoam求解器采用代數(shù)VOF(volume of fluid)方法[36]追蹤自由表面，并且在水體積分?jǐn)?shù)方程中引入了人工壓縮項(xiàng)[37]，水體積分?jǐn)?shù)輸運(yùn)方程采用MULES(Multidimensional Universal Limiter for Explicit Solution)方法進(jìn)行求解。VOF方法中引入α表示水體積分?jǐn)?shù)，α=1表示該網(wǎng)格只有水體，α=0時(shí)表示該網(wǎng)格內(nèi)只有氣體，0<α<1表示此網(wǎng)格為兩種流體的混合區(qū)域，網(wǎng)格內(nèi)存在水氣交界面，水體積分?jǐn)?shù)輸運(yùn)方程為：

(11)

在VOF水氣兩相流模型中，方程(2)中的密度和運(yùn)動(dòng)粘性系數(shù)表示為：ρ=αρwater+(1-α)ρa(bǔ)ir,μ=αμwater+(1-α)μair，ρwater、ρa(bǔ)ir分別表示水和空氣的密度，μwater、μair分別表示水和空氣的運(yùn)動(dòng)粘性系數(shù)。

2 數(shù)值驗(yàn)證與模型設(shè)置

2.1 數(shù)值驗(yàn)證

Labatut等[38]研究了MCR(Large-scale Mixed-Cell Raceway)水池中流場特性并進(jìn)行了數(shù)值模擬。本模型的有效性驗(yàn)證使用Labatut等[38]總結(jié)的MCR水池中速度分布規(guī)律，圖 1給出了在底部出流比(底面出口的流量在總流量中所占的比例)分別為0、0.15和0.20條件下速度沿著徑向的分布，在不同的底部出流比下，數(shù)值模擬結(jié)果與試驗(yàn)較為吻合；圖2給出了Labatut模擬的MCR流場和本研究使用OpenFOAM模擬的結(jié)果。根據(jù)不同工況下速度—位置分布曲線以及速度分布云圖的對比結(jié)果證明本研究的數(shù)值模型能夠有效地模擬養(yǎng)艙內(nèi)部流場。

圖1 速度沿著徑向的分布對比

圖2 MCR流場

2.2 養(yǎng)殖艙模型

圖3給出了養(yǎng)艙的幾何外形以及養(yǎng)艙在養(yǎng)殖工船中的布置。其中養(yǎng)艙艙長22.4 m、寬19.6 m、高19 m，養(yǎng)艙的底面距船體基線2 m，工作液位相對于船體基線高17 m，對應(yīng)艙內(nèi)水位為15 m。養(yǎng)艙帶有進(jìn)出水循環(huán)系統(tǒng)，進(jìn)出水口分布如圖4所示。本文選取其中一個(gè)養(yǎng)艙進(jìn)行研究，定義此養(yǎng)艙出水口的坐標(biāo)為(0,0,2)，x向?yàn)榇^方向，y向?yàn)闄M向(以左舷方向?yàn)檎?，z向?yàn)榇瓜颉?/p>

圖3 養(yǎng)艙外形以及在養(yǎng)殖工船中的位置示意圖

圖4 養(yǎng)艙進(jìn)出水口示意圖

2.3 計(jì)算設(shè)置與邊界條件

艙體壁面使用無滑移邊界條件，進(jìn)水口為速度入口，出口為壓力出口。時(shí)間步進(jìn)采用歐拉隱式，對流項(xiàng)的離散采用linearUpwind格式，壓力速度耦合求解采用pimple算法。

本研究所涉及的艙內(nèi)水的流動(dòng)是典型的內(nèi)流問題，其雷諾數(shù)：

Re=ρuL/μ

(12)

式中：ρ、μ分別為20℃時(shí)水的密度和運(yùn)動(dòng)粘性系數(shù)分別為998.2 kg/m3、1.01×10-3kg/(m·s)；u、L為流動(dòng)的特征速度和養(yǎng)艙特征尺寸，量級分別為1 m/s和10 m，可以得到此流動(dòng)的雷諾數(shù)在107量級，因此本研究使用湍流模型，采用的湍流模型為k-ωSSTDES模型。

網(wǎng)格采用六面體網(wǎng)格，在進(jìn)水口和出水口處加密，在頂部氣相區(qū)域使用較為稀疏的網(wǎng)格來降低整體網(wǎng)格量，圖 5給出了網(wǎng)格的劃分方式。

圖5 計(jì)算網(wǎng)格

2.4 網(wǎng)格收斂性

設(shè)定進(jìn)水流量為5 304 m3/h，本研究首先對網(wǎng)格收斂性進(jìn)行驗(yàn)證。網(wǎng)格采用局部加密，3種網(wǎng)格的最小分辨率分別為24 cm、12 cm和6 cm，網(wǎng)格數(shù)量分別在100萬、400萬、1 500萬左右。圖6給出了z=15.4 m的水平面上(最高的進(jìn)水口所在的平面)速度沿著x軸的分布，粗糙的網(wǎng)格的結(jié)果和精細(xì)的網(wǎng)格的差距較大，而中等網(wǎng)格的結(jié)果和比較密的網(wǎng)格的結(jié)果比較接近，最大速度出現(xiàn)的位置幾乎一致，速度最大值達(dá)到精細(xì)網(wǎng)格的96%。綜合考慮精度和算力，計(jì)算采用中等粗細(xì)的網(wǎng)格。

圖6 網(wǎng)格收斂性驗(yàn)證

3 計(jì)算結(jié)果分析

3.1自由表面渦流場特性

圖7 自由面現(xiàn)象

圖8 速度分布圖

圖9給出了本數(shù)值模擬的結(jié)果與Odgaard理論速度分布進(jìn)行對比，結(jié)果表明本研究的數(shù)值方法可以有效地捕捉到自由表面漩渦流場的速度分布。

圖9 數(shù)值模擬的速度分布與Odgaard理論對比

圖8給出了艙內(nèi)不同水平面的速度矢量圖以及x=0位置處的速度云圖示意圖，圖 10分別給出了x=0 m剖面的速度分布(漩渦流場的總速度，水平速度和垂向速度)和不同水平面上的速度分布，圖 11給了z=15.4 m平面上的速度云圖。從圖 8～圖11中可以得到在養(yǎng)艙中得到了漩渦流場，除中部的主漩渦外，在邊角處未有明顯的漩渦出現(xiàn)，可以認(rèn)為流場中沒有死水區(qū)的產(chǎn)生。和各種漩渦理論一樣，不同水平面上的水平速度幾乎一致；垂向速度在豎直方向上會(huì)隨著距離出水口的距離的增大而減小，而在遠(yuǎn)離中心軸線位置，豎直速度分布和空間位置的相關(guān)性較小。

注：(a)總速度; (b)水平速度; (c)豎直速度

圖11 z=15.4m的平面上速度分布

圖12給出了Q準(zhǔn)則識別出的Q=0.5渦量圖，在自由表面渦的軸線附近速度較大，中軸線附近的渦比較集中；渦主要分布在養(yǎng)艙中部軸線和壁面附近，遠(yuǎn)離中心漩渦和養(yǎng)艙壁面位置處渦強(qiáng)度較小。

圖12 Q=0.5的渦量圖

3.2 進(jìn)水口對流場特性影響研究

3.2.1 進(jìn)水口角度對流場特性的影響

為探究水流流入角度對速度的影響規(guī)律，本研究使用了5種不同的進(jìn)水口角度，其進(jìn)水口方向和參數(shù)如圖 13和表1所示，其中表 1 中的角度度數(shù)表示進(jìn)水口水流射入方向與x方向的夾角，不同的入射角度保證流量同為5 304 m3/h，5種工況的入射方向依次遠(yuǎn)離中心出水口，其中工況1的入射方向朝向養(yǎng)艙出水口。

表1 不同工況的4個(gè)入水口的角度(角度制)

圖13 不同工況下進(jìn)口射流方向

圖14給出了產(chǎn)生漩渦的工況5和未產(chǎn)生漩渦的工況1之間的速度分布對比，可以看到兩種工況流場存在顯著差異，與工況5的速度分布不同，工況1的速度空間分布比較均勻，速度值在一個(gè)比較小的范圍，除了出口處附近的速度與工況5較為接近外，其余大部分空間的速度和工況5的差別較大。

注：左圖為工況5；右圖為工況1；(a)x=0 m剖面水平速度; (b)x=0 m剖面垂向速度；(c)z=10 m平面水平速度

圖15(a)、15(b)分別給出10 m平面上不同工況速度沿x軸的分布，水平速度和總速度的大小接近。漩渦的產(chǎn)生對流場整體的速度提升有著積極的意義，且速度分布和進(jìn)水口的角度存在相關(guān)性，當(dāng)入水口偏離中心的角度增大時(shí)其整體的速度也會(huì)變大。漩渦流場中的速度大小會(huì)隨著進(jìn)水口的角度的改變而改變，但是一旦形成漩渦，其渦核半徑rm的大小保持穩(wěn)定。圖 15(c)為中軸線上的垂向速度沿著豎直方向的分布，垂向速度的分布在一定程度上可以反映出不同水層之間的換水效率，從圖中可以得到從自由表面到出水口，漩渦流場中垂向速度從上往下逐漸增大，表明水在一層一層自上而下匯聚，且豎直速度的大小和進(jìn)水口的角度還存在著相關(guān)關(guān)系，而未產(chǎn)生漩渦的工況1的豎直速度幾乎為0且保持不變。在出水口附近，垂向速度變化劇烈，此時(shí)未產(chǎn)生漩渦的工況1的豎直速度較大，而漩渦流場的豎直速度較小，且隨著進(jìn)水口方向偏離中心的程度的增加。出水口軸線上垂向速度減小，表明漩渦的產(chǎn)生對排水效率會(huì)產(chǎn)生一定的影響，在水利工程中漩渦會(huì)使排水效率下降是避免漩渦產(chǎn)生的一個(gè)重要原因[39]，但是在漁業(yè)養(yǎng)殖中，可以在控制同樣流量的情況下，水位相對更高，給水產(chǎn)品更充足的生存空間。

注：(a)10 m平面速度沿x軸分布；(b)10 m平面水平速度沿x軸分布；(c)垂向速度沿z軸分布

表2給出了流場的一些特征值，在同等流量下，艙內(nèi)水體積會(huì)隨著入水口的方向偏離中心程度的增大而增大。

圖17給出了速度不同進(jìn)水口角度下的概率分布，最大概率速度會(huì)隨著入水口方向偏離中心而增大，同時(shí)速度概率分布的均勻度也會(huì)隨著入水口方向偏離中心幅度的增大而變得更加均勻。為表征整體速度的大小，引入速度的概率分布函數(shù)特征值，F(xiàn)(v)為速度v的概率分布函數(shù)，v80為速度概率分布函數(shù)F(v)=0.80時(shí)對應(yīng)的速度值，計(jì)v99為速度概率分布函數(shù)F(v)=0.99對應(yīng)的速度值。

圖17 不同角度工況下速度的概率分布以及累積概率密度

表2給出了不同工況下流場的漩渦特征以及速度統(tǒng)計(jì)特征量，v80和v99均會(huì)隨著進(jìn)水口角度偏移中心的幅度的增大而增大；且隨著進(jìn)水口角度偏離中心的幅度的增大，表面凹陷的深度增大。

表2 不同角度工況下的漩渦場特性

數(shù)值計(jì)算結(jié)果表明，不同的角度工況中除工況1外，其余4種工況均產(chǎn)生中心漩渦，漩渦的強(qiáng)度會(huì)隨著進(jìn)水口方向偏離出水口程度的增大而增強(qiáng)。5種工況僅僅是進(jìn)水口角度不同，由此可知漩渦的產(chǎn)生與進(jìn)水口流入水體的推流作用有強(qiáng)烈的相關(guān)性?？紤]到漩渦形成的因素，根據(jù)角動(dòng)量守恒[40]，進(jìn)水口流入的水質(zhì)點(diǎn)在向出口處匯聚的時(shí)候切向速度會(huì)增加，進(jìn)水口方向的朝向偏離出水口的程度越大，流入的水體的初始的角動(dòng)量也就越大，其對進(jìn)水口附近的水體形成的推流作用也更加強(qiáng)烈，當(dāng)水匯聚到出口附近時(shí)的速度也會(huì)比較大。單位時(shí)間內(nèi)流入水體的角動(dòng)量和最大切向速度的關(guān)系如圖16所示，漩渦流場的速度和進(jìn)水口流入流體的角動(dòng)量正相關(guān)。

圖16 進(jìn)水口單位時(shí)間角動(dòng)量和漩渦流場最大切向速度的關(guān)系

3.2.2 進(jìn)水流量對流場特性的影響

為量化進(jìn)水流量對流場特性的影響，采用工況5的進(jìn)水方向，研究了不同流量下流場速度特性，表3給出了5種工況的流量、艙內(nèi)水位。

表3 不同流量工況下的水位

在表3所給的工況下，均都有漩渦產(chǎn)生，而且漩渦流場的速度概率分布和流量密切相關(guān)，圖18給出不同流量下速度統(tǒng)計(jì)規(guī)律。流量較大的工況流場的優(yōu)勢流速大，速度概率分布更加均勻，流速的分布范圍更廣，而流量越小其速度分布越是集中在一個(gè)較小的范圍內(nèi)。

圖18 不同流量下流速的概率分布以及累積概率密度

4 結(jié)論

本研究使用OpenFOAM建立了數(shù)值養(yǎng)艙，研究了“國信1號”養(yǎng)艙進(jìn)出水口推流影響下養(yǎng)艙流場特性。數(shù)值養(yǎng)艙復(fù)現(xiàn)了養(yǎng)艙中出現(xiàn)的漩渦流，分析了漩渦流場的速度分布特性，發(fā)現(xiàn)漩渦的產(chǎn)生對整個(gè)流場的速度有著較大的提升。研究了不同的進(jìn)水口入射角度和流量對流場的影響，數(shù)值結(jié)果顯示，隨著進(jìn)水口入射角度偏離中心程度的增大，流場的漩渦強(qiáng)度增強(qiáng)，流速增加，速度的概率分布更加均勻；固定進(jìn)水口入射角度的情況下，通過改變流量發(fā)現(xiàn)，隨著進(jìn)水口流量的增大，艙內(nèi)流場流速增加，速度的概率分布更加均勻。改變流量與進(jìn)水口入射角度都是通過改變單位時(shí)間內(nèi)流入水體的角動(dòng)量來影響流場，可以通過調(diào)整進(jìn)水口角度與進(jìn)水流量兩種方式來調(diào)整養(yǎng)艙內(nèi)的整體流速來控制適漁流速。根據(jù)該研究結(jié)果，對水體新鮮度較敏感的水產(chǎn)品養(yǎng)殖可以采用改變?nèi)肷浣嵌鹊姆绞絹斫档驼w流速，而對水體要求相對較低的水產(chǎn)品養(yǎng)殖中，可以考慮減小流量降低整體流速。