水平管降膜蒸發(fā)器殼程流場和溫度場數(shù)值研究

侯昊，畢勤成，張曉蘭

(西安交通大學 動力工程多相流國家重點實驗室，陜西 西安，710049)

淡水資源匱乏已成為當今世界威脅人類生存的嚴重問題。在地球的水資源總量中海水占97.5 %，從某種意義上說，海水是取之不盡、用之不竭的主要水源。因而，發(fā)展海水淡化技術，向海洋要淡水是解決淡水資源供需矛盾的有效途徑。海水淡化也就成為開發(fā)淡水資源的重要途徑之一。目前，已商業(yè)化的海水淡化方法有多種[1-2]，其中，低溫多效蒸餾海水淡化技術不僅熱能利用率高、產(chǎn)水率高，而且節(jié)能，近年來發(fā)展迅速，裝置規(guī)模日益擴大。另外，由于采用廉價材料，成本也日益降低，使其成為較為理想的海水淡化方式。因此，這種海水淡化技術備受各國學者的關注，成為當前海水淡化領域的研究熱點之一[3-5]。作為低溫多效蒸餾海水淡化系統(tǒng)的重要裝置，水平管降膜蒸發(fā)器受管內(nèi)蒸汽和管外海水的熱工參數(shù)影響，殼程流場處于復雜的三維流動狀態(tài)，每一局部區(qū)域都經(jīng)歷變質(zhì)量的流動、熱力和傳熱過程，并且流動與換熱相互耦合，因此整個熱力過程錯綜復雜。另外，幾何結(jié)構復雜，又封閉在殼體內(nèi)部，無法觀察和測量蒸發(fā)器內(nèi)部所有位置的流場數(shù)據(jù)，一方面使研究者很難全面了解其內(nèi)部工作過程，準確掌握整體以及局部區(qū)域真實流動和傳熱狀況的細觀信息；另一方面，又制約了試驗數(shù)據(jù)對于蒸發(fā)器結(jié)構改進和運行參數(shù)優(yōu)化的指導作用。隨著多物理場耦合技術的發(fā)展，分布參數(shù)法因其能夠描述實際物理過程，獲得更加真實準確的參數(shù)分布，受到越來越多的重視。國內(nèi)外一些學者已經(jīng)將該方法應用于某些換熱器的熱力性能分析和優(yōu)化設計[6-10]，而將其用于水平管降膜蒸發(fā)器的研究卻鮮有報道。為此，本文作者構建了水平管降膜蒸發(fā)器實體的三維分布參數(shù)模型，對蒸發(fā)器殼程區(qū)域的流場和溫度場特性進行了數(shù)值研究，分析和總結(jié)了該區(qū)域內(nèi)流體流動和傳熱的細觀規(guī)律，直觀地刻畫了蒸發(fā)器內(nèi)部的熱力過程，為深入了解復雜流動與傳熱現(xiàn)象提供了理論指導。

1 物理模型

圖1所示為水平管降膜蒸發(fā)器示意圖，其中管束分2個管程布置，1管程布置在下部，2管程在上部。換熱過程分為管內(nèi)蒸汽冷凝和管外海水蒸發(fā)2部分：一方面，加熱蒸汽首先進入1管程并大部分冷凝，剩余蒸汽經(jīng)管箱折返后進入2管程直至完全冷凝；另一方面，海水通過噴頭分布到管束上，自上而下降膜流動并逐漸吸收熱量，由過冷變?yōu)轱柡蜖顟B(tài)后開始蒸發(fā)。產(chǎn)生的二次蒸汽由管束周邊逸出，然后穿過絲網(wǎng)除沫器；在去除蒸汽中夾帶的海水液滴后，進入殼體兩側(cè)的蒸汽通道，并沿管束軸向流出。分管程布管在熱力學上是合理的，可以減小冷熱流體傳熱溫差，減少熵增引起的熱損失，另外，管內(nèi)蒸汽壓力高于管外二次蒸汽壓力，即使傳熱管因腐蝕穿透也能保證產(chǎn)品水水質(zhì)。

圖1 水平管降膜蒸發(fā)器示意圖Fig.1 Schematic of horizontal-tube falling-film evaporator

2 數(shù)學模型

2.1 基本假設

為了使模型推導方便，模型建立基于以下假設：(1) 忽略管排間液體的飛濺，即上排管底部落下的液量等于下排管頂部的初始液量；(2) 汽液界面處于熱力學平衡狀態(tài)，不考慮蒸汽的剪切力作用；(3) 不考慮傳熱管的污垢熱阻。

2.2 計算單元劃分

對于蒸發(fā)器內(nèi)部不同區(qū)域的管束而言，管內(nèi)蒸汽干度沿管子軸向逐漸減小，管外噴淋海水流量沿管排方向變化，從而造成蒸發(fā)器的換熱性能沿管子軸向和管排方向均有變化。管束示意圖如圖2所示。將蒸發(fā)器內(nèi)管束從i，j，k3個方向分別劃分為Nx，Ny，Nz等份，其中i，j，k分別表示水平、垂直和軸向方向，這樣整個管束劃分為Nx×Ny×Nz個計算單元。這種方法的優(yōu)點在于數(shù)學模型的尺度基于物理模型的實際拓撲結(jié)構，具有明確的物理意義，避免陷入CFD方法過于繁瑣的計算，極大地提高計算速度和加快程序收斂。針對每一單元建立傳熱與流動控制方程，使該單元出口參數(shù)成為下一單元的入口參數(shù)，采用迭代法求解出所有單元進出口的溫度、壓力等熱力參數(shù)。每一計算單元的熱負荷為：

式中：hoveral，hi和ho分別為總傳熱系數(shù)、管內(nèi)和管外傳熱系數(shù)；A為傳熱面積；λ為導熱系數(shù)；Ti和To分別為管內(nèi)、外工質(zhì)熱力學溫度；di和do分別為管內(nèi)、外直徑。

圖2 管束示意圖Fig.2 Schematic of tube bundle

為了計算每一單元的熱負荷和確定出口干度，需考慮蒸汽沿管子軸向和海水沿管排方向各熱工參數(shù)的變化，將物性、傳熱、流動等條件實現(xiàn)局部化，從而統(tǒng)一管內(nèi)外側(cè)相變與非相變情況下的計算，分別獲得管內(nèi)、外傳熱系數(shù)分布。另外，管內(nèi)蒸汽的壓力、溫度變化與兩相壓降有關，因而還需計算管內(nèi)汽液兩相壓降。

2.3 管內(nèi)、外傳熱系數(shù)

2.3.1 管內(nèi)傳熱系數(shù)

蒸汽進入管束時為飽和狀態(tài)，隨著熱交換的進行，經(jīng)歷汽液兩相流和單相液體狀態(tài)。當蒸汽處于汽液兩相時，Chato[11]在Nusselt冷凝換熱理論基礎上進行解析，得到了水平圓管周向平均傳熱系數(shù)，Jaster等[12]又對其進行了修正，由此得到管內(nèi)傳熱系數(shù)公式如下：

式中：ρL和ρG分別為工質(zhì)液相和氣相的密度；r為汽化潛熱；μL為液相動力黏度；g為重力加速度；Tw為管子內(nèi)壁面熱力學溫度；x為汽液兩相流干度。

當管內(nèi)蒸汽完全冷凝后，管內(nèi)傳熱系數(shù)采用齊德-泰勒公式[13]進行計算：

式中：l為管子長度；μw的定性溫度為管子內(nèi)壁面熱力學溫度。

2.3.2 管外傳熱系數(shù)

海水進入蒸發(fā)器為過冷狀態(tài)，沿管排降膜并逐漸吸熱后達到飽和狀態(tài)。對于管外液膜處于過冷狀態(tài)，采用如下公式[14]進行計算：

對于管外液膜處于飽和狀態(tài)，采用如下公式[15]進行計算：

2.4 兩相壓降

水平管內(nèi)汽液兩相流動壓降由重位壓降 Δpstatic、加速壓降Δpmom和摩擦壓降Δpfrict組成，即：

由于是水平管，重位壓降Δpstatic=0，加速壓降為：

式中：Gtotal為液汽兩相流總質(zhì)量流速；ε為空泡份額；下標in和out表示計算單元的進口和出口。

摩擦壓降采用Friedel公式[16]進行計算：

式中：ΔpL為液相摩擦壓降；φ2為兩相摩擦倍率。

3 程序計算方法

根據(jù)上述計算模型，再加上對熱力過程其他一些現(xiàn)象的數(shù)學描述及求解，開發(fā)了水平管降膜蒸發(fā)器三維穩(wěn)態(tài)過程的數(shù)值計算程序，能夠模擬不同幾何結(jié)構和流體物性參數(shù)條件下蒸發(fā)器內(nèi)部流動與傳熱過程，具體計算步驟如下。

步驟 1 輸入已知條件：包括蒸發(fā)器幾何結(jié)構參數(shù)(傳熱管規(guī)格以及管排布置結(jié)構)、加熱蒸汽和供入海水各熱工參數(shù)和相關物性參數(shù)；

步驟2 對蒸發(fā)器管束劃分計算單元；

步驟3 預估管外飽和壓力、2管程蒸汽進口壓力和流量；

步驟4 依次求解2管程和1管程各排管的熱工參數(shù)，針對任一排管，先按照軸向順序分別計算各個單元，然后沿管排方向轉(zhuǎn)至下一排管進行計算，直到蒸發(fā)器的底排管；

步驟5 根據(jù)算得的管外飽和壓力、2管程蒸汽進口壓力和流量代替舊值，重新迭代計算直至收斂。

4 計算結(jié)果及分析

4.1 模型驗證

為了驗證模型及程序的正確性與可靠性，本文對某低溫多效海水淡化系統(tǒng)的水平管降膜蒸發(fā)器進行了模擬，計算中一共劃分單元26 970個，滿足網(wǎng)格獨立性要求，并將計算結(jié)果與該蒸發(fā)器的實際運行工況進行比較。圖3所示為計算值與實際值的對比。由圖3可見：在加熱蒸汽和海水參數(shù)相同的條件下，本文計算所得二次蒸汽參數(shù)與實際工況吻合良好，兩者之間的相對誤差在-4.2%～1.6 %范圍內(nèi)，證明了本文所用模型具有很高的準確度。

4.2 海水流場及噴頭布置的影響

圖3 計算值與實際值的對比Fig.3 Comparison between calculated and actual results

圖4 海水流場分布Fig.4 Distribution of seawater flow field

圖5 改進噴頭布置后海水流場分布Fig.5 Distribution of seawater flow field for changing placement of spray nozzles

圖4和圖5所示為不同噴頭布置方式下蒸發(fā)器內(nèi)海水流場云圖，直觀地給出了海水噴淋密度沿管排方向分布的完整數(shù)據(jù)。由圖4和5可知：蒸發(fā)器內(nèi)噴頭布液效果良好，降膜傳熱過程順暢進行。由圖4可知：蒸發(fā)器采用均勻噴頭布置方式時，海水噴淋密度由上部的63 g/(ms)減小到下部的 40 g/(ms)，呈現(xiàn)出“上高下低”的分布規(guī)律，其原因在于過冷的海水需要2管程蒸汽加熱后才能達到飽和狀態(tài)，而且由于2管程蒸汽量較少，其換熱量也較少，即使海水已經(jīng)飽和，產(chǎn)汽量也非常小，因此海水噴淋密度在2管程時大體上保持不變。到達1管程后，海水開始蒸發(fā)，產(chǎn)生二次蒸汽。在1管程中部海水噴淋密度由60 g/(ms)逐漸降低到46 g/(ms)，而且整體變化較為平穩(wěn)，說明此區(qū)域換熱效果良好，產(chǎn)汽量保持穩(wěn)定；而在1管程進口的下半?yún)^(qū)海水噴淋密度已經(jīng)很小，其最小值為40 g/(ms)，但仍能夠保證傳熱管束周向完全潤濕，不會引起液膜破裂。同時，在該區(qū)域海水噴淋密度沿管長方向由40 g/(ms)逐漸增大到50 g/(ms)，表現(xiàn)出分布不均勻的特點，不利于避免“干斑”現(xiàn)象發(fā)生。因此，可以采用“前密后疏”的噴頭布置方式來改善這一區(qū)域的海水流動情況。如圖5所示，2管程上部區(qū)域海水分布不均，噴淋密度對應地由 70 g/(ms)逐漸減小到 60 g/(ms)，而1管程下部區(qū)域海水分布則更加均勻，噴淋密度大體上保持不變，約為48 g/(ms)。

4.3 海水溫度場

圖6所示為蒸發(fā)器內(nèi)海水溫度場云圖，過冷的海水由50 ℃逐漸達到飽和溫度57 ℃，說明2管程已經(jīng)達到了對海水的預熱作用。在2管程出口的上半?yún)^(qū)海水溫度始終保持在50 ℃，這是因為管內(nèi)蒸汽已經(jīng)完全冷凝，兩相換熱變?yōu)閱蜗鄵Q熱，過冷海水無法得到足夠熱量使其溫度明顯上升，還要繼續(xù)吸熱才能達到飽和狀態(tài)。根據(jù)能量守恒原理，加熱蒸汽帶入的熱量等于進出口海水變化的熱量與二次蒸汽帶走的熱量之和，因此，可以通過對海水進行預熱提高海水的初始溫度來有效減小其過冷度，進而增大二次蒸汽所占的熱量份額，提高加熱蒸汽的熱利用率，從而增加單效蒸發(fā)器的二次蒸汽產(chǎn)量。

圖6 海水溫度場分布Fig.6 Distribution of seawater temperature field

4.4 傳熱特征

蒸發(fā)器內(nèi)部的傳熱特征如圖7～9所示。由圖7可知：蒸汽進入管束后沿管長方向逐漸冷凝，通過釋放潛熱將熱量傳給管外海水，干度也逐漸降低且變化較為平穩(wěn)，說明蒸汽沿管長流動過程中并未出現(xiàn)急劇冷凝現(xiàn)象，整體傳熱狀況良好。在1管程出口處平均干度約為0.15，剩余蒸汽經(jīng)管箱折返后進入2管程直至完全冷凝。另外，總傳熱系數(shù)也可以表征各個計算單元的換熱效果，與管內(nèi)、外傳熱系數(shù)密切相關，特別受蒸汽質(zhì)量流速、干度以及海水噴淋密度影響較大。由圖8和圖9可知：在1管程和2管程蒸汽進口處總傳熱系數(shù)最大，約為3.3 kW/(m2K)，而后沿管長方向逐漸減小，說明進口區(qū)換熱較強，特別是2管程表現(xiàn)得更加明顯，因為該區(qū)域蒸汽干度大，管內(nèi)擾動強烈，使得管內(nèi)傳熱系數(shù)較大；而后沿管長方向冷凝換熱逐漸減弱，總傳熱系數(shù)也發(fā)生較大變化，管內(nèi)單相流動時總傳熱系數(shù)很小，僅相當于管內(nèi)蒸汽冷凝情況下總傳熱系數(shù)的7%左右。在1管程進口的下半?yún)^(qū)總傳熱系數(shù)略有減小，原因在于上半?yún)^(qū)海水已經(jīng)達到飽和狀態(tài)，開始蒸發(fā)汽化，噴淋密度沿管排方向逐漸減小，管外液膜波動減緩，削弱了管外換熱，進而影響了總體傳熱。此外，整個蒸發(fā)器總傳熱系數(shù)的平均值為2.8 kW/(m2K)，已經(jīng)達到了設計參數(shù)要求。

圖7 蒸汽干度分布Fig.7 Distributions of steam quality

圖8 2管程總傳熱系數(shù)分布Fig.8 Distributions of overall heat transfer coefficient for second tube pass

圖9 1管程總傳熱系數(shù)分布Fig.9 Distributions of overall heat transfer coefficient for first tube pass

4.5 海水鹽度分布

蒸發(fā)器內(nèi)海水鹽度的變化情況如圖10所示。與圖4所表述的流場分布一一對應，表現(xiàn)出“上小下大”的趨勢。由于2管程海水的噴淋密度基本保持不變，其鹽度也維持在36 g/kg。當海水到達1管程后開始蒸發(fā)，產(chǎn)生二次蒸汽，噴淋密度逐漸減小，鹽度則隨之升高。在1管程中部海水鹽度整體變化較為平穩(wěn)，從36 g/kg逐漸增加到48 g/kg。在1管程進口的下半?yún)^(qū)海水鹽度明顯升高，最大值達到56 g/kg，再次說明此區(qū)域換熱強烈，生成較多的二次蒸汽。但是海水鹽度過高，會促使其在蒸發(fā)過程中所含鹽分受熱分解并轉(zhuǎn)化成不同形態(tài)物質(zhì)沉淀在傳熱管表面形成結(jié)垢，因此在工程設計中應當格外引起注意，控制海水濃縮比不超過2.5，最高溫度低于70 ℃可以避免發(fā)生上述現(xiàn)象。

圖10 海水鹽度分布Fig.10 Distributions of seawater salinity

5 結(jié)論

(1) 建立了水平管降膜蒸發(fā)器實體的三維分布參數(shù)模型，數(shù)值解與實際值的相對誤差在-4.2%～1.6%范圍內(nèi)，驗證了該模型的可靠性，該模型適用于單元復疊循環(huán)計算且成本低廉。

(2) 通過模擬蒸發(fā)器內(nèi)部的熱力過程，獲得了大量流場、溫度場和傳熱參數(shù)分布信息，真實地反映了蒸發(fā)器內(nèi)部復雜的流動與傳熱特性，以細觀角度揭示出局部區(qū)域流體流動和傳熱的深層次問題。

(3) 三維分布參數(shù)模型是進行水平管降膜蒸發(fā)器數(shù)值研究的一種合理且實用的模型，對蒸發(fā)器的機理研究及整體的優(yōu)化設計具有指導意義。

[1]周赤忠,李焱. 當前海水淡化主流技術的分析與比較[J]. 電站輔機,2008,29(4) : 1-5.ZHOU Chi-zhong,LI Yan. Analysis and comparison between current major technologies of seawater desalination[J]. Power Station Auxiliary Equipment,2008,29(4): 1-5.

[2]宿翠霞,倪華秋,李凌霄,等. 海水淡化技術及其應用現(xiàn)狀[J].中國資源綜合利用,2009,27(7) : 31-32.SU Cui-xia,NI Hua-qiu,LI Ling-xiao,et al. Seawater desalination technology and its application[J]. China Resources Comprehensive Utilization,2009,27(7): 31-32.

[3]Alasfour F N,Darwish M A,Bin Amer A O. Thermal analysis of ME-TVC+MEE desalination systems[J]. Desalination,2005,174(1): 39-61.

[4]Darwish M A,El-Dessouky H. The heat recovery thermal vapor-compression desalting system: A comparison with other thermal desalination processes[J]. Applied Thermal Engineering,1996,16(6): 523-537.

[5]Dey P K,Hawlader M N A,Chou S K,et al. Performance of a single-effect desalination system operating with different tube profiles and materials[J]. Desalination,2004,166(15): 69-78.

[6]Zhang L N,Yang C X,Zhou J H. A distributed parameter model and its application in optimizing the plate-fin heat exchanger based on the minimum entropy generation[J]. International Journal of Thermal Sciences,2010,49(8): 1427-1436.

[7]Rodriguez M,Diaz M S. Dynamic modeling and optimisation of cryogenic systems[J]. Applied Thermal Engineering,2007,27(7):1182-1190.

[8]Wang F Q,Maidment J F,Missenden J F,et al. A novel special distributed method for dynamic refrigeration system simulation[J]. International Journal of Refrigeration,2007,30(5):887-903.

[9]Sultana P,Wijeysundera N E,Ho J C,et al. Modeling of horizontal tube-bundle absorbers of absorption cooling system[J].International Journal of Refrigeration,2007,30(4): 709-723.

[10]Zavala-Rio A,Santiesteban-Cos R. Reliable compartmental models for double-pipe heat exchangers: An analytical study[J].Applied Mathematical Modelling,2007,31(9): 1739-1752.

[11]Chato J C. Laminar condensation inside horizontal and inclined tubes[J]. American Society of Heating Refrigerating and Airconditioning Engineers Journal,1962,4: 52-60.

[12]Jaster H,Kosky P G. Condensation heat transfer in a mixed flow regime[J]. International Journal of Heat and Mass Transfer,1976,19(1): 95-99.

[13]楊世銘,陶文銓. 傳熱學[M]. 北京: 高等教育出版社,1998:168-169.YANG Shi-ming,TAO Wen-quan. Heat transfer[M]. Beijing:Higher Education Press,1998: 168-169.

[14]Sernas V. Heat transfer correlation for subcooled water films on horizontal tubes[J]. Transactions of the ASME,1979,101(1):176-178.

[15]Parken W H,Fletcher L S,Sernas V,et al. Heat transfer through falling-film evaporation and boiling on horizontal tubes[J].Journal of Heat Transfer,1990,112(3): 744-750.

[16]Ould Didi M B,Kattan N,Thome J R. Prediction of two-phase pressure gradients of refrigerants in horizontal tubes[J].International Journal of Refrigeration,2002,25(7): 935-947.