漁獲物冰漿保鮮流動傳熱分析

趙新穎，黃溫赟，黃文超，呂續(xù)艦

(1 中國水產(chǎn)科學研究院漁業(yè)機械儀器研究所，農(nóng)業(yè)農(nóng)村部遠洋漁船與裝備重點實驗室，上海 200092;2 海洋國家實驗室深藍漁業(yè)工程裝備技術(shù)聯(lián)合實驗室，山東 青島 266237;3 南京理工大學能源與動力工程學院，江蘇 南京 210094)

冰漿是一種由平均直徑約0.10 mm的冰粒與水溶液組成的混合物[1]，具有良好的流動性和較大的蓄冷密度，是目前廣泛應用的蓄冷二次循環(huán)制冷介質(zhì)之一[2-3]。為降低建造成本，國內(nèi)漁船大部分未安裝壓縮式制冷裝置，一般攜帶碎冰給漁獲物保鮮，效果不佳，影響漁獲物品質(zhì)。利用冰漿的流動傳熱對漁獲物進行降溫冷藏，確保漁獲物品質(zhì)安全，是一種切實可行的方法。研究冰漿的流動和傳熱特性，掌握船載冰漿循環(huán)理論，具有重要現(xiàn)實意義。

國內(nèi)學者在冰漿領域開展了大量頗具成效的研究工作。江煥寶[4]通過對某極地運輸船海水系統(tǒng)中作為滑油冷卻器的管殼式換熱器的數(shù)值模擬，實現(xiàn)了換熱器的結(jié)構(gòu)參數(shù)優(yōu)化。劉圣春等[5]總結(jié)了換熱器的流動換熱特性研究現(xiàn)狀，認為確定兩相顆粒流動的流變模型是解決目前換熱問題的前提。文力等[6]對冰漿制取裝置的過冷器中流體流動進行了模擬分析。楊帆等[7]利用Bingham模型分析了冰漿在層流和湍流狀態(tài)下的摩擦因子,得到了冰漿在水平直管內(nèi)流動的壓降曲線。熊庭等[8]對水平管道中泥漿輸送的固液兩相流進行了CFD數(shù)值模擬，并與Durand 模型的計算值和文獻試驗數(shù)據(jù)進行了對比分析。徐愛祥等[9]通過建立 CFD-PBM 耦合模型，研究含冰率與流速對水平直管內(nèi)冰漿的冰晶體積分數(shù)與粒徑的分布規(guī)律。鄢新[10]采用Fluent對冰水固流轉(zhuǎn)化的傳熱過程進行了仿真模擬，根據(jù)模擬結(jié)果初步給出了冰水界面熱傳導問題的解決辦法。葉健[11]通過試驗研究了質(zhì)分數(shù)小于19.2%的TBAB水合物漿體在內(nèi)徑為2.0 mm和4.5 mm圓直管中的流動特性，通過測量水合物漿體的流速和壓降得到了相應的流動曲線。

此外，何國庚等[12]早在2005年就注意到了冰漿的應用前景，并對冰漿流體的應用研究提出了一些有價值的建議。王繼紅等[13-14]采用計算流體力學方法研究了水平管道內(nèi)冰漿流體的無相變流動過程，給出了管道內(nèi)冰漿流體阻力特性；鄧義斌等[15]基于商用CFD軟件Fluent對水平管內(nèi)冰水兩相流動進行了計算研究，發(fā)現(xiàn)管內(nèi)兩相流動阻力特性受顆粒碰撞彈型恢復系數(shù)、水流流速、顆粒直徑、顆粒質(zhì)量流量等特性的影響較為顯著。張鵬等[16]同樣考慮水平管內(nèi)的冰漿流動情況，不同之處在于考慮了冰漿在三相點處溫度下水平管內(nèi)的流動情況；白銀等[17]通過試驗對不同時間、不同加熱量、不同流速下傳熱器進出口的冰漿的溫度、含冰率進行了測試；徐立等[18-19]研究了冰漿在極地船上的應用，包括極地船換熱器中海水-冰晶兩相流的流動及傳熱特性，以及冰晶顆粒在極地船殼管式換熱器海水管內(nèi)的分布和融化特性[20-21]。

相較于碎冰等傳統(tǒng)方式，采用冰漿對漁獲物進行保鮮能夠提供更為均勻的冷量分布，但相關研究工作鮮有報道。本研究建立了冰漿流動的計算模型，首先采用文獻數(shù)據(jù)對計算模型進行了驗證，隨后研究了入口冰顆粒體積分數(shù)對漁獲物上表面與側(cè)面的熱通量特性的影響，最后研究了流動速度對漁獲物表面?zhèn)鳠崃康挠绊?，力圖為冰漿在漁獲物保鮮中的應用提供具有一定實用價值的參考。

1 數(shù)學模型

1.1 模型概述

采用歐拉-歐拉雙流體模型對冰漿流動傳熱進行仿真計算，其基本思想是將液體相和顆粒相等效為可以互相貫穿的連續(xù)性介質(zhì)，針對液體相和顆粒擬流體相分別建立N-S方程，同時液體相和顆粒相之間的耦合是通過相間作用力來實現(xiàn)的。顆粒間的相互作用以碰撞為主，壁面處通過設置鏡面反射系數(shù)和彈性恢復系數(shù)來表現(xiàn)顆粒與壁面的相互作用。其連續(xù)性方程為[16]：

(1)

式中：α、ρ和v分別為兩相的體積分數(shù)、密度(kg/m3)及速度矢量(m/s)；i表示l(液相)或者s(固相)，下文中的式(3)和(5)～(7)中下標q表示與i相對的相(i=l，q=s；i=s，q=l)。

1.2 動量守恒方程

動量守恒方程可以表示為[16]：

(2)

式中：τ為應力張量(N/m2)；P為流場壓力分布；g是重力加速度(m/s2)；FD,i和FL,i分別表示冰顆粒在第i相中的拖曳力和升力(N)；為哈密頓算子。對于液相有：

(3)

而對于固相，根據(jù)擬流體假設，其應力張量τS由顆粒碰撞產(chǎn)生的隨機粒子運動引起：

(4)

式中：μS、ξS和PS分別表示剪切黏度(N·s/m2)、體積黏度(N·s/m2)以及固相壓強(N·s/m2)；為單位向量。

1.3 相間作用力

在歐拉雙流體模型中，兩相間的耦合關系是通過相間作用力建立起來的，本次計算所考慮的相間作用力包括拖曳力和升力。拖曳力是顆粒在液相中的阻力，是固液兩相間最基本的相互作用形式[16]，其表達式為:

FD,i=Ksl(vq-vi)

(5)

式中：Ksl表示顆粒的拖曳系數(shù)，可采用Gidaspow模型進行計算：

(6)

顆粒在液相中運動時還受到由液相速度梯度引起的升力，其表達式為：

FL,i=CLαSρl(|vl-vs|)×(×vl)

(7)

式中：CD和CL分別為顆粒阻力系數(shù)和升力系數(shù)。

結(jié)合上述模型，考慮顆粒動力學理論、標準k-ε湍流模型和能量方程模型，即可實現(xiàn)冰漿流動和傳熱特性分析[16]。

2 水平管內(nèi)冰漿流動模型驗證

2.1 初始條件及邊界條件

初始時水平管內(nèi)冰顆粒均勻分布，速度垂直于入口截面，且固相與液相速度相同。管道出口采用壓力出口條件，出口為一個大氣壓。壁面采用絕熱條件，液相無滑移，冰顆粒在壁面處采用Johnson-Jackson模型處理，顆粒與壁面碰撞恢復系數(shù)為0.9，鏡面系數(shù)為0.015。

2.2 網(wǎng)格劃分及計算方法

首先在ANSYS-ICEM中建立了水平管道三維模型，管道直徑D為23 mm，為了保證流動充分發(fā)展，管道長度L設置為2 000 mm。管道網(wǎng)格采取了六面體結(jié)構(gòu)網(wǎng)格，網(wǎng)格數(shù)為72萬個，采用非穩(wěn)態(tài)求解模型，求解方法為壓力與速度耦合的SIMPLE求解器，控制方程采用QUICK與二階迎風離散格式。開始計算時時間步長設置為較小的0.000 5，以保證收斂性；每個時間步長內(nèi)控制方程的殘差小于10-4時認為計算收斂。

2.3 計算結(jié)果及驗證

在Fluent軟件中對水平管內(nèi)冰漿兩相流動進行數(shù)值求解，計算了在冰顆粒體積分數(shù)為0.1、速度為1.0 m/s時的冰顆粒體積分數(shù)分布情況，并與Zhang等[1]計算結(jié)果進行驗證，以保證利用這種方法模擬冰漿與漁獲物流動傳熱的可行性與可靠性，如圖1所示，吻合情況較好。由于冰顆粒的密度略小于海水，因此在流動過程中冰顆粒會不斷向上懸浮，使得出口截面上的冰顆粒體積分數(shù)自上而下遞減，圓管頂部最高為0.103 5，底部最低為0.09。

圖1 水平管出口截面冰顆粒體積分數(shù)驗證(U=1.0 m/s)

3 方管內(nèi)冰漿與漁獲物的流動傳熱計算

3.1 物理模型及計算條件

在冰漿的管道輸送中，外界的傳熱容易使冰漿在傳輸過程中發(fā)生相變。根據(jù)漁獲物尺寸，在ANSYS ICEM中建立方管模型并劃分網(wǎng)格。如圖2所示，方管截面為邊長4 m的正方形，長度20 m；管內(nèi)漁獲物的長、寬、高分別為3.2 m、2.36 m、2.1 m，網(wǎng)格數(shù)為48萬個。

圖2 方管內(nèi)冰漿和漁獲物流動傳熱模型

由圖2可知，漁獲物在方管內(nèi)沿x軸向的位置為x=10～13.2 m，漁獲物頂部在截面y=0.1上，兩側(cè)面分別位于z=-1.18和z=1.18截面上。

采用歐拉-歐拉雙流體模型對模型進行計算，在初始時刻方管內(nèi)冰顆粒均勻分布，固相與液相溫度均為273.15 K；冰漿速度垂直于入口截面，且固相與液相速度相同；管道出口設置為壓力出口，條件為一個大氣壓；管內(nèi)漁獲物邊界條件設置為wall，給定壁溫298.15 K，且壁面處考慮冰顆粒的碰撞與反射。采用耦合SIMPLE求解器進行穩(wěn)態(tài)計算，控制方程采用二階迎風方法進行離散。

3.2 計算結(jié)果

3.2.1 不同冰顆粒體積分數(shù)對傳熱效果的影響

熱通量(Heat flux)又稱為熱流，表示單位時間通過單位面積的熱能，可以用傅里葉定律來充分描述。與溫度相關的熱通量可以描述為：

(8)

式中：k為表面熱導率；τ為溫度分布；x為沿表面的某一坐標，負號表示熱通量從高溫區(qū)傳向低溫區(qū)。在對漁獲物與流體傳熱進行計算時，利用了冰顆粒與漁獲物表面之間的熱通量來表征冰漿對漁獲物的冷卻情況，在CFD-Post 19.0中讀取漁獲物表面?zhèn)鳠崃?，全管如圖3所示。

圖3 熱通量分布情況(U=0.5 m/s，αin=0.2)

冰顆粒的分布是影響傳熱的主要原因之一，圖4展示了流速為0.5 m/s，冰顆粒入口體積分數(shù)αin=0.2時漁獲物附近截面上冰顆粒分布情況。

圖4 漁獲物附近冰顆粒分布情況

從圖4中可以看出，在冰漿流動充分發(fā)展以后，由于冰的密度較小，在重力作用下冰顆粒大量聚集在方管頂部，在方管出口截面(x=20)處冰顆粒體積分數(shù)最高為0.397；管道中間部分則較為均勻，體積分數(shù)在0.2～0.22；而在漁獲物附近，方管底部幾乎沒有冰顆粒分布；在冰漿繞流漁獲物時，由于渦的產(chǎn)生以及壁面處速度梯度較大，湍流脈動較強，導致漁獲物附近冰顆粒質(zhì)量濃度較低，且分布不均勻，繞流結(jié)束后冰顆粒分布有重新趨于穩(wěn)定的趨勢。

冰漿質(zhì)量濃度主要由冰顆粒體積分數(shù)Vf體現(xiàn)，傳熱效果與冰顆粒的分布情況密切相關。當一定流速的冰漿在方管中沖擊漁獲物迎流面后會發(fā)生繞流，冰顆粒與漁獲物表面之間發(fā)生碰撞反射，與顆粒自身的彈性碰撞以及冰水兩相相互作用，必然導致冰顆粒分布不均勻，使得漁獲物各表面?zhèn)鳠嵝Ч煌?，當冰漿繞流漁獲物的頂部與側(cè)面時，這種情況尤為明顯。故分別計算了αin=0.1、0.2、0.3時，漁獲物上表面(x=10)與側(cè)面(z=1.18)冰顆粒體積分數(shù)Vf與熱通量[q(W/m2)]的分布情況，并讀取上表面和側(cè)面中線上的數(shù)據(jù)，其變化趨勢曲線如圖5所示。

圖5 不同αin下漁獲物上表面與側(cè)面熱通量分布情況

圖6中漁獲物上表面冰顆粒質(zhì)量濃度沿z軸的變化較為明顯，其總體呈中心低、兩側(cè)高的規(guī)律分布，且隨著入口冰顆粒體積分數(shù)αin的提高，表面熱通量呈現(xiàn)出整體增大的趨勢。圖7展示出漁獲物表面附近冰漿的流線分布情況?？梢钥闯觯?/p>

圖6 不同αin情況下漁獲物上表面和側(cè)面熱通量和冰顆粒體積分數(shù)變化趨勢

圖7 漁獲物上表面附近流線分布

(1)上表面中心處熱通量較小，原因可能為當冰漿繞流90°棱角的漁獲物時，上表面位于方管中心位置，流速較高，顆粒與壁面碰撞所受的排斥力較大，使得湍流脈動較為強烈，冰漿繞流后在漁獲物表面上方形成渦流，使冰顆粒有遠離壁面的趨勢，進而影響了顆粒與壁面間的傳熱，而在αin=0.3的工況下，發(fā)現(xiàn)冰顆粒體積分數(shù)與表面熱通量的波動明顯小于另外兩種工況；

(2)在上表面兩側(cè)冰顆粒的質(zhì)量濃度和熱通量都相對較高，αin為0.1、0.2和0.3時熱通量最高分別達到9 221.9 W/m、21 396.5 W/m和21 375.3 W/m。因為在冰漿繞流過程中，由于冰顆粒密度較低，且繞過棱角時速度高，壓力較低，造成漁獲物兩側(cè)的流體向上表面翻轉(zhuǎn)的現(xiàn)象，因此上表面兩側(cè)冰顆粒質(zhì)量濃度較高，熱通量明顯大于中心位置；

(3)這種向上翻轉(zhuǎn)的情況同樣使得漁獲物側(cè)面的熱通量分布不均，側(cè)面與上表面接觸的位置熱通量明顯高于其他位置，在αin為0.2時最高達到21 612.3W/m，超過αin為0.1時的2倍；

(4)在重力作用下，冰顆粒隨著流動過程逐漸向管道頂部聚集，使得漁獲物側(cè)面的冰顆粒體積分數(shù)與熱通量自上而下呈遞減趨勢。

3.2.2 不同流速對傳熱效果的影響

對于固液兩相流動而言，在近壁面區(qū)域，冰顆粒所受到的升力和顆粒與壁面之間碰撞所產(chǎn)生的排斥力對其分布情況影響較大；而升力與液相速度梯度成正比，顆粒與壁面的碰撞作用主要受到入口體積分數(shù)αin的影響，因此在相同入口體積分數(shù)下，排斥力隨著入口速度的增大而增大，故管道上壁面聚集的冰顆粒隨著入口速度的增大呈減小的趨勢，冰顆粒會逐漸趨向于管道中間位置，如圖8所示。對αin=0.2，入口速度分別為0.5 m/s、1.5 m/s和2.5 m/s的情況下漁獲物表面冰顆粒分布情況以及傳熱情況計算發(fā)現(xiàn)：

圖8 不同流速下漁獲物表面冰顆粒體積分數(shù)與熱通量變化曲線

(1)隨著流動速度的增大，漁獲物表面冰顆粒體積分數(shù)和熱通量均有所提高，其中U=1.5 m/s情況下，在靠近側(cè)棱約0.2 m的距離內(nèi)熱通量略高于2.5 m/s的工況，且側(cè)棱處冰顆粒體積分數(shù)基本相同，其原因可能為在速度較高時，側(cè)棱處發(fā)生繞流時冰顆粒慣性較大，導致和漁獲物表面距離較大，造成熱通量反轉(zhuǎn)的情況；

(2)在流速為1.5 m/s和2.5 m/s的工況下，漁獲物側(cè)面熱通量在靠近側(cè)棱處的冰顆粒質(zhì)量濃度和熱通量相差無幾；相較于上表面，漁獲物側(cè)面冰顆粒體積分數(shù)和熱通量的變化較大，在2.5 m/s工況下漁獲物側(cè)面最高熱通量約是0.5 m/s工況下的4.5倍，分別為104 822 W/m和23 132 W/m，原因可能是在流速較低的情況下冰顆粒易于向方管頂部聚集，導致繞流漁獲物的下層冰漿冰顆粒體積分數(shù)較小，因此傳熱效果相對較差。

4 結(jié)論

冰顆粒的分布情況是影響其與漁獲物表面?zhèn)鳠岬闹饕蛩?，隨著冰顆粒入口體積分數(shù)的提高，漁獲物上表面與側(cè)面的熱通量整體有所增大；隨著入口冰顆粒體積分數(shù)的增大，其在漁獲物表面的分布趨于均勻。通過對3種流動速度的工況進行對比可以發(fā)現(xiàn)，適當提高流動速度可以增大漁獲物表面冰顆粒的體積分數(shù)，從而增加表面熱量交換。當流動速度增大時，漁獲物表面的冰顆粒數(shù)量一開始明顯增大，當速度增大到一定水平時，漁獲物表面上冰顆粒體積分數(shù)十分接近，此時漁獲物上表面?zhèn)鳠崃吭龇苄　?/p>

□