侯昊，畢勤成，張曉蘭

(西安交通大學(xué) 動(dòng)力工程多相流國(guó)家重點(diǎn)實(shí)驗(yàn)室，陜西 西安，710049)

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

1 物理模型

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

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

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

2.1 基本假設(shè)

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

2.2 計(jì)算單元?jiǎng)澐?/h3>

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

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

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

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

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

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

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

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

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

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

2.3.2 管外傳熱系數(shù)

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

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

2.4 兩相壓降

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

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

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

摩擦壓降采用Friedel公式[16]進(jìn)行計(jì)算：

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

3 程序計(jì)算方法

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

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

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

步驟3 預(yù)估管外飽和壓力、2管程蒸汽進(jìn)口壓力和流量；

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

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

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

4.1 模型驗(yàn)證

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

4.2 海水流場(chǎng)及噴頭布置的影響

圖3 計(jì)算值與實(shí)際值的對(duì)比Fig.3 Comparison between calculated and actual results

圖4 海水流場(chǎng)分布Fig.4 Distribution of seawater flow field

圖5 改進(jìn)噴頭布置后海水流場(chǎng)分布Fig.5 Distribution of seawater flow field for changing placement of spray nozzles

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

4.3 海水溫度場(chǎng)

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

圖6 海水溫度場(chǎng)分布Fig.6 Distribution of seawater temperature field

4.4 傳熱特征

蒸發(fā)器內(nèi)部的傳熱特征如圖7～9所示。由圖7可知：蒸汽進(jìn)入管束后沿管長(zhǎng)方向逐漸冷凝，通過釋放潛熱將熱量傳給管外海水，干度也逐漸降低且變化較為平穩(wěn)，說明蒸汽沿管長(zhǎng)流動(dòng)過程中并未出現(xiàn)急劇冷凝現(xiàn)象，整體傳熱狀況良好。在1管程出口處平均干度約為0.15，剩余蒸汽經(jīng)管箱折返后進(jìn)入2管程直至完全冷凝。另外，總傳熱系數(shù)也可以表征各個(gè)計(jì)算單元的換熱效果，與管內(nèi)、外傳熱系數(shù)密切相關(guān)，特別受蒸汽質(zhì)量流速、干度以及海水噴淋密度影響較大。由圖8和圖9可知：在1管程和2管程蒸汽進(jìn)口處總傳熱系數(shù)最大，約為3.3 kW/(m2K)，而后沿管長(zhǎng)方向逐漸減小，說明進(jìn)口區(qū)換熱較強(qiáng)，特別是2管程表現(xiàn)得更加明顯，因?yàn)樵搮^(qū)域蒸汽干度大，管內(nèi)擾動(dòng)強(qiáng)烈，使得管內(nèi)傳熱系數(shù)較大；而后沿管長(zhǎng)方向冷凝換熱逐漸減弱，總傳熱系數(shù)也發(fā)生較大變化，管內(nèi)單相流動(dòng)時(shí)總傳熱系數(shù)很小，僅相當(dāng)于管內(nèi)蒸汽冷凝情況下總傳熱系數(shù)的7%左右。在1管程進(jìn)口的下半?yún)^(qū)總傳熱系數(shù)略有減小，原因在于上半?yún)^(qū)海水已經(jīng)達(dá)到飽和狀態(tài)，開始蒸發(fā)汽化，噴淋密度沿管排方向逐漸減小，管外液膜波動(dòng)減緩，削弱了管外換熱，進(jìn)而影響了總體傳熱。此外，整個(gè)蒸發(fā)器總傳熱系數(shù)的平均值為2.8 kW/(m2K)，已經(jīng)達(dá)到了設(shè)計(jì)參數(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所表述的流場(chǎng)分布一一對(duì)應(yīng)，表現(xiàn)出“上小下大”的趨勢(shì)。由于2管程海水的噴淋密度基本保持不變，其鹽度也維持在36 g/kg。當(dāng)海水到達(dá)1管程后開始蒸發(fā)，產(chǎn)生二次蒸汽，噴淋密度逐漸減小，鹽度則隨之升高。在1管程中部海水鹽度整體變化較為平穩(wěn)，從36 g/kg逐漸增加到48 g/kg。在1管程進(jìn)口的下半?yún)^(qū)海水鹽度明顯升高，最大值達(dá)到56 g/kg，再次說明此區(qū)域換熱強(qiáng)烈，生成較多的二次蒸汽。但是海水鹽度過高，會(huì)促使其在蒸發(fā)過程中所含鹽分受熱分解并轉(zhuǎn)化成不同形態(tài)物質(zhì)沉淀在傳熱管表面形成結(jié)垢，因此在工程設(shè)計(jì)中應(yīng)當(dāng)格外引起注意，控制海水濃縮比不超過2.5，最高溫度低于70 ℃可以避免發(fā)生上述現(xiàn)象。

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

5 結(jié)論

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

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

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

[1]周赤忠,李焱. 當(dāng)前海水淡化主流技術(shù)的分析與比較[J]. 電站輔機(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]宿翠霞,倪華秋,李凌霄,等. 海水淡化技術(shù)及其應(yīng)用現(xiàn)狀[J].中國(guó)資源綜合利用,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]楊世銘,陶文銓. 傳熱學(xué)[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.