齊春華，徐　克，馮厚軍，呂慶春

(1.國家海洋局天津海水淡化與綜合利用研究所，天津 300192；2.國家海水利用工程技術研究中心，天津 300192)

低溫多效蒸餾海水淡化主要采用水平管降膜蒸發(fā)技術，水平管降膜蒸發(fā)器傳熱傳質性能直接影響著海水淡化成本。長期以來，作為熱交換的核心元件，傳熱管均是以圓形截面為主，隨著計算流體力學和微觀傳熱技術的發(fā)展，人們越來越認識到橢圓管對液膜流動和強化傳熱的影響，并成功應用于空調與制冷、化工等行業(yè)，表現(xiàn)出較好的節(jié)能效果。國內外許多專家學者對橢圓管換熱開展了多方面的研究[1-4]，2003年，Hasan[5]對橢圓管和圓管進行了試驗研究，表明橢圓管的傳熱與流阻性能是圓管的1.93～1.96倍。在換熱設計中，殼程阻力降是很重要的，它直接影響工藝流程上下游的穩(wěn)定性、泵功耗和操作費用[6]。在低溫多效海水淡化過程中，管外液膜分布狀態(tài)將直接影響傳熱效率以及噴淋供液泵的耗能。采用橢圓管換熱元件，對改變管外薄膜流動狀態(tài)、提升降膜蒸發(fā)效率有積極作用，但橢圓管降膜過程中的管外液膜流動相關研究較少，其在低溫多效海水淡化中的應用研究也未見報道。

本研究采用開發(fā)后的FLUENT對橢圓管外液膜流動進行了三維數(shù)值模擬研究，得出管型對液膜分布、液膜厚度的影響，為開發(fā)基于橢圓管的降膜蒸發(fā)海水淡化技術做探索性研究。

1　機理分析

水平管降膜蒸發(fā)是20世紀70年代發(fā)展起來的一項蒸餾淡化技術，液體從分布器噴淋到管子頂部，即成薄膜狀沿管子兩側流下進行降膜蒸發(fā)，直到管底匯聚成液滴滴到下一排管子上，如圖1所示。因為傳熱管兩側均有相變發(fā)生，換熱系數(shù)較高，而管外壁降膜流動狀態(tài)是蒸發(fā)側傳熱膜系數(shù)的決定因素。

圖1　圓管噴淋降膜物理模型Fig.1　Laminar flow physical-model of spraying liquid films on horizontal tube

圖1為水平管外噴淋降膜的物理模型，在忽略二次蒸汽對液膜的剪切力作用時，可以認為液體流動的推動力是重力沿管壁的切向分力Gg，Gg隨著周向角θ增大有明顯的變化。在0<θ<π/2區(qū)域內，Gg隨著θ增大逐漸增大，液膜運動速度逐漸加快，造成上游液體供應不上，液膜減薄；在θ=π/2處，Gg等于重力G，液體受到推動力最大，流速達到最大值，液膜厚度出現(xiàn)最小值；在π/2<θ<π區(qū)域內，Gg隨著θ增大逐漸減小，液膜的運動速度也不斷減小，由于下游液體流速比上游慢，故在下游出現(xiàn)液體堆積，液膜厚度不斷增加，直至θ=π處，Gg=0，液體脫離管壁向下滴?？梢?，在θ=π/2處的附近區(qū)域，液膜流速大，膜厚小，擾動加劇，故熱阻較小，傳熱系數(shù)較高。若使這個較高傳熱系數(shù)的區(qū)域所占比例變大，傳熱效果會得到強化，所以截面形狀為橢圓型的傳熱管就是根據(jù)上述基礎理論開發(fā)的。在形狀因子E(長短軸之比)大于1的橢圓管周長中的較大區(qū)域(簡稱為B區(qū))上，作用在液膜上的推動力較大，而且這個區(qū)域隨形狀因子E的增大而增大。圖2是Moalen等[7]通過數(shù)學推導得出的結果。在相同液膜雷諾數(shù)和傳熱溫差條件下，E=2.0的橢圓管B區(qū)占半周長的75%，而E=0.5的橢圓管和E=1.0的圓管分別只占14%和54%；當E>2.0時，B區(qū)范圍很大，噴淋蒸發(fā)傳熱膜系數(shù)he很大，而E<2.0(E=0.5和1.0)的傳熱管噴淋蒸發(fā)傳熱性能要差得多。

圖2　傳熱膜系數(shù)的周向分布[7]Fig.2　Variation of the heat transfer coefficient in circumferential direction[7]

2　數(shù)值模擬方法

采用計算流體力學方法對橢圓管型建立液膜流動模型進行數(shù)值模擬，確定管外液膜流動狀態(tài)較好地管型結構尺寸。

2.1　物理模型建立

水平管降膜蒸發(fā)器橢圓管單管外流動換熱模型如圖3所示，傳熱管上方有小孔布液器。傳熱管頂部正對布液器開孔中心，過布液器開孔中心的縱剖面作為研究截面。利用scheme語言技術進行參數(shù)化設計；在Gambit中針對橢圓管型分別進行參數(shù)化建模、網(wǎng)格劃分和邊界條件設置；在Fluent中進行了參數(shù)化初始條件設置、邊界條件設置、預處理、求解和后處理。

圖3　橢圓管單管外流動換熱模型Fig.3　Flow and heat transfer model outside the single elliptical tube

2.2　基本假設及流體熱物性

假設所選取的計算區(qū)域內流體的流動為湍流，計算流體為自來水，處于蒸發(fā)器內壓力為0.027 MPa的工作環(huán)境，假設蒸發(fā)器內部除了計算流體外充滿了空氣，流體是不存在換熱情況下的純流動，其物性參數(shù)為常數(shù)，如表1所示。

表1　流體物性參數(shù)Table 1　Thermophysical parameters of liquid

2.3　控制方程和VOF方法

傳熱管外汽液兩相流動遵循的基本方程組，包括動量、質量和能量守恒方程和雷諾方程，以及為封閉雷諾方程提出的各種湍流模型等。通過對一定區(qū)域內基本方程組進行分析，就可以得到在特定情況下流體流動的規(guī)律。其通用形式：

(1)

Φ=1，ΓΦ=0，SΦ=0

(2)

得出連續(xù)性方程。

(3)

得出動量方程。

(4)

得出雷諾方程。

采用模擬兩相流通用的VOF方法來追蹤氣液相界面。在VOF表達式中，動量方程與前面的表達式相同，而且求解方法也沒有變化。但連續(xù)方程需要稍作修正。對于某一相有如下表達式：

?α2/?t+v·α2=0

(5)

式(5)中，α2為第2相得體積分數(shù)。第1相的體積分數(shù)為：α1=1-α2。流體性質的估算就可以通過體積分數(shù)來表達，例如，在任何網(wǎng)格單元的物質密度可用式(6)來表達：

ρ=α2ρ2+(1-α2)ρ1

(6)

式(6)中，ρ是在某一計算單元格內的估算密度，ρ1和ρ2分別為第1相和第2相的密度。

2.4　網(wǎng)格劃分和邊界條件

為了提高計算精度，對計算區(qū)域網(wǎng)格的處理采用分區(qū)劃分、局部加密的方式。整個計算域長(x)×高(y)×寬(z)=18 mm × 50 mm × 40 mm，分成壁面區(qū)、臨界壁面區(qū)、入流區(qū)、軸向區(qū)、出流區(qū)和非主要區(qū)，共6個區(qū)，全部采用六面體網(wǎng)格，網(wǎng)格總數(shù)為1 420 000個，計算區(qū)域網(wǎng)格如圖4所示。

圖4　網(wǎng)格劃分圖Fig.4　Meshed configuration of the bath

建立相應的邊界條件如下：1)計算區(qū)域邊線包括2種邊界：速度入口(Velocity inlet)邊界及Wall邊界，速度入口邊界位于計算區(qū)域上邊線，寬度為1.5 mm； 2)Wall邊界為傳熱管外壁，其余部分為對稱邊界。模擬過程中計算區(qū)域初始給定為全氣相，布液孔入流液體為全液相。

3　模擬結果與分析

3.1　管型對液膜分布的影響

在入口流速分別為0.2、0.5和0.8 m/s[液體負荷分別為0.037、0.093和0.149 kg/(m·s)]的工況下，對相同截面周長的圓管型EP01和橢圓管型EP03進行管外液體流動的數(shù)值模擬計算。

表2　EP01和EP03管型參數(shù)Table 2　The structural parameters of the EP01 and EP03 tubes

圖5為圓管型不同入口流速下的液膜分布情況，其中圖5a)為EP01圓管型不同入口流速下的液膜分布情況，圖5b)為EP03橢圓管型不同入口流速下的液膜分布情況。

圖5　不同入口流速下的液膜分布情況Fig.5　Falling film distribution chart with the velocity

從圖5a)可見，液滴在每個管型頂部完成滴落、鋪展、下落和脫離再到恢復初始狀態(tài)的全過程，其脫落時間隨著入口流速的增大而減少。從圖5b)可見，橢圓管上的液膜在流速為0.2 m/s時，液體呈不連續(xù)的滴狀模式流動，傳熱管表面存在大量干壁區(qū)域；流速為0.5 m/s時，管外液體流型呈柱狀模式，但液體不能完全潤濕傳熱管外壁，仍存在干壁區(qū)域；流速為0.8 m/s時，管外液體流型為滴-柱狀模式，一層薄液膜包裹在傳熱管外壁，而當流速增加至1.0 m/s時，管外液體流型呈現(xiàn)出柱狀模式，管外包裹了一層較厚的液膜，并且2個入流口中間區(qū)域的液膜厚度明顯增加，形成明顯的隆起。

由圖5a)和圖5b)對比可見，液膜在圓管底端脫落的分布較為零亂，液膜在管底部聚集，脫落較慢，使管底部液膜厚度大大增加，從圖5b)入口流速0.8 m/s時可見，液膜在橢圓管底端脫落呈間隔柱狀，均勻程度比圓管有較大程度的改善。橢圓管外液膜局部流速增大，薄膜效果增強，且很快從管底部脫落，基本沒有液膜聚集現(xiàn)象。

由圖5可見，橢圓管外液膜的形成、鋪展、下落和脫離所用時間都比圓管短，入口流速0.8 m/s時，橢圓管用時335 ms，而圓管用時375 ms，表明橢圓管外液體流速快，更新快，擾動增強，利于換熱。

3.2　管型對液膜厚度的影響

液體在管表面的分布決定著管外液膜厚度，而傳熱管管型對管外液體的分布有重要影響。在液體負荷為0.149 kg/(m·s)，即入口流速0.8 m/s時，根據(jù)模擬結果對管外壁液膜厚度沿管周和管軸方向的變化情況進行了測量，測量結果見圖6，其中圖6a)為管外液膜厚度分布，圖6b)為管長方向液膜厚度分布(θ=90°)，θ為測點和管重心的連線與過管重心的垂線之間的夾角。

由圖6a)可見，液膜在傳熱管表面的開始分布及流動階段，在大部分區(qū)域，圓管的膜厚比橢圓管厚10%～20%，而只有在管周液體流速最大處(θ=90°)，圓管外液膜厚度才與橢圓管外液膜厚度基本相當；液膜快脫落時，在管底端的液膜厚度，圓管達3.6～3.8 mm，而橢圓管也達2.7～2.9 mm，分別是θ=90 °時的11和8倍；橢圓管僅在最低點厚度較厚，在其他區(qū)域則沒有任何液膜聚集現(xiàn)象，橢圓管液膜厚度在整個管周區(qū)域均低于圓管。

由圖6b)為在過管軸的平面與管周相交處(θ=90 °)液膜厚度沿管長方向的變化，可見橢圓管外的液膜較圓管外薄，在L=9 mm處膜厚最厚，這是因為2個入流口(0 mm處和18 mm處)中間區(qū)域會因液膜的相向鋪展而造成液體堆積，但由圖6可見橢圓管外的液體局部堆積現(xiàn)象明顯減輕。

圖6　液膜厚度分布Fig.6　Film thickness vs.variation

可見，橢圓管相對于圓管，其管外壁液膜脫落速度加快，分布均勻，易形成穩(wěn)定柱狀流，管表面的液膜薄，更利于傳熱。

3.2.1同一流速下不同橢圓管管外液膜厚度對比

在入口流速為0.8 m/s的工況下，對物理模型就截面周長相等而截面形狀因子不同的幾種橢圓管進行管外液體流動的數(shù)值模擬計算，結構尺寸如下表：

表3　5種不同橢圓管管型參數(shù)Table 3　The structural parameters of five kind of elliptical tubes

由圖7可見，截面形狀因子E越大，橢圓管外膜厚越小，但隨著E的不斷增大，液膜的不均勻性也越來越嚴重(見圖8)。如當E=2.1時，θ=170 °和θ=20 °處，橢圓管外液膜厚度分別為0.55 mm、0.16 mm，相差2倍多。E=1.5要比E=1.2的液膜薄且更加均勻。

圖7　不同橢圓管外的液膜厚度變化Fig.7　The falling film thickness around the different elliptical tube

圖8　EP05橢圓管外液膜厚度分布圖Fig.8　Vaper-liquid phase distribution chart of the EP05 elliptical tube

通過模擬研究，E=1.5的橢圓管具有液膜薄，均勻性好的特點，管外液體流動形態(tài)更易形成柱狀流，有利于傳熱。

3.2.2同一橢圓管不同入口流速橢圓管外液膜厚度對比

對E=1.5的橢圓管，在入口流速為0.2～1.4 m/s(θ=90 °)時進行管外液體流動的數(shù)值模擬計算。

圖9　入口流速的變化對橢圓管外側液膜水平處厚度影響Fig.9　Variation of falling film thickness around elliptical tube with the velocity

由圖9可見，隨著入口流速的增大，橢圓管外的液體負荷也增加，管外液膜的厚度相應增大。入口流速為0.2 m/s時，模擬結果顯示管外存在較大面積的干壁，降低蒸發(fā)效率；當入口流速大于0.8 m/s時，液膜厚度增加1.5倍左右，影響傳熱，且提高供應液體泵的功率消耗，不利于節(jié)能。

模擬結果表明，橢圓管形狀因子E=1.5，入口流速V=0.8 m/s時，管外液膜狀態(tài)最好，液膜薄且均勻，利于換熱。

4　結論

1)首先建立橢圓管外降膜流動的物理模型，通過模擬研究管外液體流動狀態(tài)，發(fā)現(xiàn)橢圓管型相對于圓管型，液體流動耗費時間縮短40 ms，說明橢圓管型上的表面液體更新速度比圓管加快，利于換熱。

2)在傳熱管表面大部分區(qū)域的液膜，圓管的膜厚比橢圓管厚10%～20%；在管長方向上，橢圓管外的液膜厚度較圓管外薄，且局部聚集現(xiàn)象明顯減輕；截面形狀因子E越大，橢圓管外膜厚越小，但液膜的不均勻性也越嚴重。

3)模擬結果表明，E=1.5的橢圓管具有液膜薄，均勻性好的特點，管表面更易形成柱狀流，當入口流速V=0.8 m/s時，此橢圓管外液膜狀態(tài)最好。

參考文獻：

[1]劉希女，楊永安，鄒同華，等.叉排橢圓管強迫對流換熱的數(shù)值模擬[J].制冷與空調，2004(1):16-19

[2]LI Q，MA L，QI S，etal. Experimental investigation on heat transfer and f low resistance of an 1.6/1.0 elliptic cylinder in cross-flow[J].Journal of Qingdao University of Science and Technology，2004，25(1)：61-54

[3]LI Q，QI S，ZHANG X，etal.Heat transfer characteristics and flow behaviors of elliptic cylinders in cross-flow[J].Journal of Qingdao University of Science and Technology，2004，25(5)：434-440

[4]SABOYA S，SABOYA F E M.Experiments on elliptic sections in one-and tow-row arrangements of plate fin and tube heat exchangers [J].Experimental Thermal and Fluid Science，2001，24：67-75

[5]HASAN A，SIREN K.Performance investigation of plain and finned tube evaporatively cooled heat exchangers [J].Applied Thermal Engineering，2003，23 (3)：325-340

[6]HORVAT A，LESKOVAR M，MAVKO B.Comparison of heat transfer conditions in tube bundle cross-flow for different tube shapes[J].International Journal of Heat and Mass Transfer，2006，49(6)：1 027-1 038

[7]MOALEM D，SIDEMAN S.Theoretical analysis of a horizontal condenser-evaporator tube[J].International Journal of Heat and Mass Transfer，1976，19(3)：259-270

[8]齊春華，李炎,俞永江,等.水平管降膜蒸發(fā)海水淡化的傳熱試驗研究[J].化學工業(yè)與工程,2011,28(5):1-5