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

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

隨著社會的發(fā)展，以豬、牛為代表的“紅肉”消費(fèi)逐漸將被以雞、魚為代表的“白肉”代替，聯(lián)合國糧農(nóng)組織(FAO)的統(tǒng)計數(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)害，被稱為“移動的海洋牧場”，工船養(yǎng)殖是發(fā)展深遠(yuǎn)海養(yǎng)殖的一個重要的方向[3-4]。

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

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

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

以上關(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ù)性方程和動量方程分別為：

(1)

(2)

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

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

(3)

(4)

其中：

(5)

(6)

(7)

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

方程(3)中湍動能k的耗散項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ù)方程中引入了人工壓縮項[37]，水體積分?jǐn)?shù)輸運(yùn)方程采用MULES(Multidimensional Universal Limiter for Explicit Solution)方法進(jìn)行求解。VOF方法中引入α表示水體積分?jǐn)?shù)，α=1表示該網(wǎng)格只有水體，α=0時表示該網(wǎng)格內(nèi)只有氣體，0<α<1表示此網(wǎng)格為兩種流體的混合區(qū)域，網(wǎng)格內(nèi)存在水氣交界面，水體積分?jǐn)?shù)輸運(yùn)方程為：

(11)

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

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

2.1 數(shù)值驗證

Labatut等[38]研究了MCR(Large-scale Mixed-Cell Raceway)水池中流場特性并進(jìn)行了數(shù)值模擬。本模型的有效性驗證使用Labatut等[38]總結(jié)的MCR水池中速度分布規(guī)律，圖 1給出了在底部出流比(底面出口的流量在總流量中所占的比例)分別為0、0.15和0.20條件下速度沿著徑向的分布，在不同的底部出流比下，數(shù)值模擬結(jié)果與試驗較為吻合；圖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所示。本文選取其中一個養(yǎng)艙進(jìn)行研究，定義此養(yǎng)艙出水口的坐標(biāo)為(0,0,2)，x向為船頭方向，y向為橫向(以左舷方向為正)，z向為垂向。

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

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

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

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

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

Re=ρuL/μ

(12)

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

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

圖5 計算網(wǎng)格

2.4 網(wǎng)格收斂性

設(shè)定進(jìn)水流量為5 304 m3/h，本研究首先對網(wǎng)格收斂性進(jì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%。綜合考慮精度和算力，計算采用中等粗細(xì)的網(wǎng)格。

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

3 計算結(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)生。和各種漩渦理論一樣，不同水平面上的水平速度幾乎一致；垂向速度在豎直方向上會隨著距離出水口的距離的增大而減小，而在遠(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個入水口的角度(角度制)

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

圖14給出了產(chǎn)生漩渦的工況5和未產(chǎn)生漩渦的工況1之間的速度分布對比，可以看到兩種工況流場存在顯著差異，與工況5的速度分布不同，工況1的速度空間分布比較均勻，速度值在一個比較小的范圍，除了出口處附近的速度與工況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)入水口偏離中心的角度增大時其整體的速度也會變大。漩渦流場中的速度大小會隨著進(jìn)水口的角度的改變而改變，但是一旦形成漩渦，其渦核半徑rm的大小保持穩(wěn)定。圖 15(c)為中軸線上的垂向速度沿著豎直方向的分布，垂向速度的分布在一定程度上可以反映出不同水層之間的換水效率，從圖中可以得到從自由表面到出水口，漩渦流場中垂向速度從上往下逐漸增大，表明水在一層一層自上而下匯聚，且豎直速度的大小和進(jìn)水口的角度還存在著相關(guān)關(guān)系，而未產(chǎn)生漩渦的工況1的豎直速度幾乎為0且保持不變。在出水口附近，垂向速度變化劇烈，此時未產(chǎn)生漩渦的工況1的豎直速度較大，而漩渦流場的豎直速度較小，且隨著進(jìn)水口方向偏離中心的程度的增加。出水口軸線上垂向速度減小，表明漩渦的產(chǎn)生對排水效率會產(chǎn)生一定的影響，在水利工程中漩渦會使排水效率下降是避免漩渦產(chǎn)生的一個重要原因[39]，但是在漁業(yè)養(yǎng)殖中，可以在控制同樣流量的情況下，水位相對更高，給水產(chǎn)品更充足的生存空間。

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

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

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

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

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

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

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

圖16 進(jìn)水口單位時間角動量和漩渦流場最大切向速度的關(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)計規(guī)律。流量較大的工況流場的優(yōu)勢流速大，速度概率分布更加均勻，流速的分布范圍更廣，而流量越小其速度分布越是集中在一個較小的范圍內(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)生對整個流場的速度有著較大的提升。研究了不同的進(jìn)水口入射角度和流量對流場的影響，數(shù)值結(jié)果顯示，隨著進(jìn)水口入射角度偏離中心程度的增大，流場的漩渦強(qiáng)度增強(qiáng)，流速增加，速度的概率分布更加均勻；固定進(jìn)水口入射角度的情況下，通過改變流量發(fā)現(xiàn)，隨著進(jìn)水口流量的增大，艙內(nèi)流場流速增加，速度的概率分布更加均勻。改變流量與進(jìn)水口入射角度都是通過改變單位時間內(nèi)流入水體的角動量來影響流場，可以通過調(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)殖中，可以考慮減小流量降低整體流速。