馮留海， 聶永廣， 王江云， 雷曉東， 王　娟， 毛　羽

(1.中國石油大學(xué) 重質(zhì)油國家重點(diǎn)實(shí)驗(yàn)室， 北京 102249； 2.北京低碳清潔能源研究所， 北京 102209；3.新奧科技發(fā)展有限公司， 河北 廊坊 065001； 4.中機(jī)康元糧油裝備(北京)有限公司， 北京 100083)

基于非平衡熱力學(xué)模型的管內(nèi)沸騰過程模擬

馮留海1,2， 聶永廣1,3， 王江云1， 雷曉東4， 王娟1， 毛羽1

(1.中國石油大學(xué) 重質(zhì)油國家重點(diǎn)實(shí)驗(yàn)室， 北京 102249； 2.北京低碳清潔能源研究所， 北京 102209；3.新奧科技發(fā)展有限公司， 河北 廊坊 065001； 4.中機(jī)康元糧油裝備(北京)有限公司， 北京 100083)

摘要：為了更準(zhǔn)確地模擬壁面沸騰過程，根據(jù)Lavievile非平衡熱力學(xué)模型對沸騰流動的壁面熱流率進(jìn)行了劃分，通過調(diào)節(jié)函數(shù)結(jié)合Sato模型修正了沸騰流動中的混合k-ε模型，建立了同時適用于管內(nèi)過冷沸騰與整體沸騰流動階段的曳力模型。使用上述模型對環(huán)狀豎直管內(nèi)沸騰流動過程進(jìn)行數(shù)值模擬，計算結(jié)果與實(shí)驗(yàn)值吻合，可以用來模擬整個管內(nèi)沸騰流動過程。數(shù)值模擬顯示，沸騰流動可以明顯增大管壁的傳熱系數(shù)，但當(dāng)管壁處蒸氣體積分?jǐn)?shù)超過0.6時，會出現(xiàn)傳熱惡化現(xiàn)象。

關(guān)鍵詞：沸騰流動； 非平衡熱力學(xué)； 曳力模型； 傳熱系數(shù)； 數(shù)值模擬

管內(nèi)沸騰流動是石油、化工、能源動力等工業(yè)領(lǐng)域的一種常見現(xiàn)象，例如石油化工行業(yè)常用的管式加熱爐和管殼式換熱器中的管內(nèi)沸騰流動。沸騰流動可以顯著提高管內(nèi)傳熱效率，所以受到人們的廣泛關(guān)注。早在20世紀(jì)60年代，Bartolomei等[1]就實(shí)驗(yàn)研究了空泡分率對垂直管內(nèi)高壓水沸騰過程的影響。隨著測量技術(shù)的發(fā)展，Bartel[2]、Lee等[3]對各種工況下的管內(nèi)沸騰過程進(jìn)行了更加詳細(xì)的測量。一方面希望實(shí)驗(yàn)獲得管內(nèi)沸騰過程中詳細(xì)的兩相流局部參數(shù)特性，另一方面希望通過實(shí)驗(yàn)測量結(jié)合理論方法歸納出氣泡脫離直徑、成核頻率和氣泡停留時間等管內(nèi)沸騰流動中的關(guān)鍵參數(shù)，從而完善數(shù)值模型，使其可以準(zhǔn)確模擬管內(nèi)沸騰過程。Kurul等[4]最早基于雙流體模型對壁面熱流率進(jìn)行了劃分，成功建立了RPI模型。因其考慮了壁面熱流率的分配情況，所以能更真實(shí)地反映管內(nèi)高壓水過冷沸騰過程，并在相關(guān)領(lǐng)域得到廣泛采用。此后，大量學(xué)者不斷完善RPI的子模型。Anglart等[5]采用RPI模型研究了核燃料分配器中的過冷沸騰流動過程，Kocar等[6-7]建立了低壓工況下水的過冷沸騰模型，李祥東等[8-9]詳細(xì)總結(jié)了面向液氮過冷沸騰流動的數(shù)值模型，竇從從等[10-11]通過數(shù)值模擬確定了垂直管內(nèi)高壓高過冷水過冷沸騰流動的子模型。研究表明，RPI模型能夠較好地預(yù)測過冷沸騰流動特性，但是由于其假設(shè)蒸氣均為飽和溫度，從而影響了管內(nèi)沸騰流動過程中整體沸騰階段的模擬精度。筆者在上述基礎(chǔ)上，并根據(jù)Lavievile的非平衡熱力學(xué)方法重新劃分了壁面熱流率[12]，通過數(shù)值計算修正了相關(guān)湍流模型和相間動量傳輸模型，最后將相關(guān)模型通過UDF植入通用CFD計算軟件Fluent，以期更好地模擬管內(nèi)沸騰流動過程。

1管內(nèi)沸騰過程的數(shù)學(xué)模型

1.1雙流體模型

雙流體模型假設(shè)氣、液兩相均為連續(xù)介質(zhì)，根據(jù)牛頓力學(xué)法則建立各自的控制方程，兩相之間的相互作用通過控制方程中的相間傳輸相描述。在管內(nèi)沸騰模型中，將液相視為連續(xù)相，蒸氣相視為擬連續(xù)相。雙流體模型是整個管內(nèi)沸騰數(shù)值模型的骨架，管內(nèi)沸騰流動數(shù)值模擬中所有子模型旨在封閉雙流體模型的控制方程組。適用于管內(nèi)沸騰流動模擬的雙流體模型的基本控制方程組為式(1)～(3)。

質(zhì)量守恒方程：

(1)

動量守恒方程：

(2)

能量守恒方程：

(3)

式中所有與相變相關(guān)的項都需要采用沸騰模型進(jìn)行計算。

1.2壁面沸騰模型

管內(nèi)沸騰是由加熱壁面加熱液體產(chǎn)生的沸騰現(xiàn)象，包括RPI模型在內(nèi)的適用于這一過程的模型總稱為壁面沸騰模型。在這一模型中，管內(nèi)的傳質(zhì)傳熱主要分成2部分，即壁面處的傳質(zhì)傳熱和主流區(qū)域的傳質(zhì)傳熱，其中計算壁面區(qū)域傳質(zhì)和傳熱的關(guān)鍵在于合理劃分壁面熱流率。

(1)壁面熱流率劃分

前人數(shù)值模擬研究采用的RPI模型大多將壁面熱流率劃分為3部分，即液相對流熱流率Ql、激冷熱流率Qq和蒸發(fā)熱流率Qe。需要說明的是，Qq是由于加熱產(chǎn)生的氣泡脫離壁面時補(bǔ)充的冷液體被瞬間加熱升溫所產(chǎn)生的熱流率。具體的表達(dá)式為式(4)～(6)。

Ql=hl(Tw-Tl)(1-Aq)

(4)

(5)

(6)

現(xiàn)有的RPI模型假設(shè)流場中的蒸氣均處于飽和溫度狀態(tài)，不考慮蒸氣在壁面處的對流換熱。本研究采用的Lavievile壁面熱流率劃分方法是在RPI模型的基礎(chǔ)上加入了蒸氣對流熱流率Qv，其表達(dá)式為式(7)。壁面熱流率的劃分和分配函數(shù)f(αl) 如式(8)和(9)所示。式(6)中的氣泡脫離直徑dbw使用Unal模型計算[13]，而主流區(qū)的dbw根據(jù)過冷度[8,10-11]確定。

Qv=hv(Tw-Tv)

(7)

Qw=f(αl)(Ql+Qq+Qe)+(1-f(αl))Qv

(8)

(9)

(2)相間傳質(zhì)和傳熱

相間質(zhì)量傳遞分為2部分。一是壁面附近液體在加熱作用下沸騰相變?yōu)檎魵猓钦魵庠谶^冷度較高的主流區(qū)域冷凝為液體。蒸氣相的傳質(zhì)項表達(dá)式為式(10)，液相與蒸氣相的相間傳熱項Qlv表達(dá)式為式(11)。

(10)

Qlv=hlv(Tl-Tv)Alv

(11)

式(10)中，hlv根據(jù)Ranz-Marshall模型計算。

與RPI模型不同，因?yàn)榭紤]蒸氣處于非飽和溫度狀態(tài)，因此既存在液相向氣相的熱量傳遞，也存在氣相向液相的熱量傳遞。在壁面附近，壁面熱流率中的對流熱流率Ql和激冷熱流率Qq最終加入到液相能量源項，蒸發(fā)熱流率Qe和蒸氣對流熱流率Qv加入到氣相能量源項。

(3)湍流模型修正

RPI模型模擬過冷沸騰過程時，通常假設(shè)液相為湍流流動，氣相為層流流動。但隨著沸騰流動過程的發(fā)展，當(dāng)管內(nèi)蒸氣體積分?jǐn)?shù)超過0.3時，這一假設(shè)明顯欠合理。因此，筆者使用混合k-ε模型計算兩相的混合湍流動能km和混合湍流耗散率εm，然后采用調(diào)節(jié)函數(shù)fα計算最終蒸氣湍流黏度。調(diào)節(jié)函數(shù)fα的表達(dá)式為式(12)。

(12)

式(12)中，a決定了函數(shù)曲線的變化幅度，b決定了函數(shù)曲線的閾值。圖1為調(diào)整函數(shù)fα曲線。從圖1可以看出，取a=20、b=0.3比較合適。最終有效黏度表達(dá)式為式(13)。

(13)

考慮到氣泡擾動對液相湍流的影響，使用Sato模型[14]修正液相有效湍流黏度，表達(dá)式為式(14)～(15)。

(14)

(15)

(4)相間曳力模型修正

在對管內(nèi)過冷沸騰流動數(shù)值模擬中，Ishii模型[15]表現(xiàn)出了良好的適用性。Ishii模型原理上僅適用于氣相體積分?jǐn)?shù)小于0.5的沸騰流動過程，而實(shí)際過程中隨著蒸氣量的增加，管內(nèi)氣液兩相流的流型則從泡狀流逐漸過渡到霧狀流。所以，有必要建立同時適用于泡狀流動和霧狀流動的曳力模型，以研究整個管內(nèi)沸騰流動過程。采用式(12)中的調(diào)節(jié)函數(shù)，將系數(shù)a、b的值分別取為20和0.5，修正的曳力系數(shù)的表達(dá)式如式(16)所示。

圖1　調(diào)整函數(shù)fα曲線

K=fα·Kbubble+(1-fα)·Kdroplet

(16)

式(16)中的Kbubble采用Ishii模型計算，而Kdroplet則采用適用性較好的Schiller模型計算。

2加熱管幾何模型及模擬工況

加熱管數(shù)值計算的幾何模型如圖2所示。模擬對象為Bartel的管內(nèi)沸騰過程[2]，其內(nèi)層為加熱棒，外層為玻璃管。加熱棒直徑為9.55 mm，加熱部分的長度為1500 mm(淺灰色部分)，玻璃管內(nèi)徑為19.05 mm。操作工況為常壓、入口溫度367.9 K，工質(zhì)為水，質(zhì)量流率為700 kg/(m2·s)，壁面熱流率為188 kW/m2，出口設(shè)置為壓力出口。因?yàn)橛嬎阌驗(yàn)檩S對稱結(jié)構(gòu)，計算在柱坐標(biāo)系中進(jìn)行，網(wǎng)格劃分采用四邊形網(wǎng)格。壓力與速度耦合采用PCSIMPLE(Phase couple SIMPLE)算法，對流項采用QUICK差分格式。

圖2　加熱管幾何結(jié)構(gòu)及尺寸

3結(jié)果與討論

圖3為Sato模型中不同Cμb取值時加熱管z=1500 mm處液相湍流黏度沿徑向的分布。由圖3看到，Cμb決定了在湍流效應(yīng)中是液體剪切造成的湍流還是氣泡誘導(dǎo)產(chǎn)生的湍流起主導(dǎo)地位。當(dāng)Cμb取值較小時，液體剪切造成的湍流占主導(dǎo)地位，這時液體湍流黏度在加熱壁面與管壁的中間位置最高，而在管壁附近則較低；壁面沸騰時管壁附近加熱產(chǎn)生的氣泡會不斷脫離壁面，破壞邊界層，增加湍流擾動，促使壁面附近的氣、液兩相快速混合，因此液相湍流黏度應(yīng)該在加熱壁面附近較高。而Cμb最適宜的取值需要進(jìn)一步分析其對管內(nèi)溫度場和蒸氣體積分?jǐn)?shù)分布造成的影響而定。

圖3　不同Cμb取值時加熱管內(nèi)z=1500 mm處液

圖4為不同Cμb時 加熱管z=1500 mm處的液體溫度的徑向分布。從圖4可以看出，Cμb取值越大，表明壁面附近液體的溫度越低，說明氣泡誘導(dǎo)的湍流作用增強(qiáng)了壁面附近的熱傳遞，從而有更多的熱量被帶到主流區(qū)域。在主流區(qū)域，由于氣相體積分?jǐn)?shù)的減小，氣泡產(chǎn)生的擾動作用減小，所以液體剪切作用產(chǎn)生的湍流重新占據(jù)主導(dǎo)地位。液體溫度會直接影響加熱壁面處氣泡脫離直徑的大小，進(jìn)而影響到蒸發(fā)熱流率的大小，因此液相湍流修正對準(zhǔn)確預(yù)測加熱管內(nèi)蒸氣體積分?jǐn)?shù)起到了非常重要的作用。

圖5為不同Cμb時加熱管z=1500 mm處的蒸氣體積分?jǐn)?shù)沿徑向的分布。從圖5可以看出，Cμb為0.6時，計算得到的徑向蒸氣體積分?jǐn)?shù)分布與實(shí)驗(yàn)值最為吻合，說明Cμb=0.6時計算得到的氣泡誘導(dǎo)湍流作用最能反映實(shí)際情況。此外，加熱管內(nèi)蒸氣體積分?jǐn)?shù)最高處并非處于管壁上，而是在臨近加熱壁面的位置。這是因?yàn)樵诒诿鏉櫥Φ淖饔孟?，加熱壁面附近生成的氣泡會被推離壁面，最終在壁面附近形成蒸氣富集區(qū)。

圖4　不同Cμb取值時加熱管內(nèi)z=1500 mm

圖5　不同Cμb時加熱管內(nèi)z=1500 mm處蒸氣體積分?jǐn)?shù)

圖6為使用修正后的曳力模型和Ishii模型計算得到的加熱管內(nèi)蒸氣體積分?jǐn)?shù)軸向分布的對比。從圖6可以看出，在過冷沸騰段，使用修正后的曳力模型計算得到的體積分?jǐn)?shù)軸向分布與Ishii模型的計算結(jié)果基本吻合，而在蒸氣體積分?jǐn)?shù)超過0.4后(過冷沸騰向整體沸騰過渡)，修正后的曳力模型計算得到的體積分?jǐn)?shù)高于Ishii模型的計算結(jié)果。修正后的曳力模型在管內(nèi)沸騰的高體積分?jǐn)?shù)區(qū)域與實(shí)驗(yàn)值吻合更好。

圖6　不同曳力模型模擬的加熱管內(nèi)蒸氣體積分?jǐn)?shù)

圖7是加熱管壁面?zhèn)鳠嵯禂?shù)與壁面附近蒸氣體積分?jǐn)?shù)軸向分布。從圖7可以看出，在加熱段開始階段，壁面?zhèn)鳠嵯禂?shù)隨著蒸氣體積分?jǐn)?shù)的增加迅速增加；在蒸氣體積分?jǐn)?shù)略大于0.1時，壁面?zhèn)鳠嵯禂?shù)達(dá)到峰值，之后存在一個迅速下降階段和一個穩(wěn)定發(fā)展階段，隨后開始呈現(xiàn)線性遞減趨勢。分析認(rèn)為，在管內(nèi)沸騰的初始階段，壁面產(chǎn)生的氣泡破壞了邊界層，增強(qiáng)了壁面附近的流體湍流作用，導(dǎo)致壁面熱傳遞系數(shù)迅速升高；當(dāng)壁面沸騰進(jìn)一步發(fā)展時，氣泡產(chǎn)生的湍流也增強(qiáng)了氣泡和液體之間的熱傳遞和質(zhì)量傳遞，而這時壁面附近液體的過冷度依然較高，使得壁面上氣泡的產(chǎn)生和壁面附近氣泡的冷凝達(dá)到短暫的平衡狀態(tài)，壁面附近蒸氣體積分?jǐn)?shù)增加趨勢放緩，同時這種平衡狀態(tài)也弱化了壁面附近的對流熱傳遞，使得壁面?zhèn)鳠嵯禂?shù)減小；管內(nèi)沸騰流動繼續(xù)發(fā)展，管壁附近液體被冷凝的氣泡充分加熱使得液體過冷度明顯降低，這時壁面附近的氣泡將不會被液體完全冷凝，加熱管壁附近的蒸氣體積分?jǐn)?shù)開始持續(xù)增加，而加熱管壁附近不斷增加的蒸氣會使得壁面?zhèn)鳠嵯禂?shù)降低。從圖7可以推斷，當(dāng)加熱管壁面附近的蒸氣體積分?jǐn)?shù)超過0.6時，壁面?zhèn)鳠嵯禂?shù)甚至?xí)陀趩蜗嘁后w加熱時的壁面?zhèn)鳠嵯禂?shù)。

圖7　加熱管壁面?zhèn)鳠嵯禂?shù)和壁面附近

4結(jié)論

(1)根據(jù)Lavievile的非平衡熱力學(xué)方法對加熱管壁面沸騰中的壁面熱流率進(jìn)行劃分，使其可以適用于整個管內(nèi)沸騰過程的數(shù)值計算。同時，對壁面沸騰中的兩相湍流模型修正，建立了同時適用于過冷沸騰和整體沸騰的雙流體相間曳力模型。數(shù)值模擬顯示，在使用混合k-ε模型時，Sato模型的修正系數(shù)Cμb取0.6時，計算結(jié)果與實(shí)驗(yàn)值吻合較好。

(2)在豎直管內(nèi)沸騰中，壁面潤滑力的作用使得蒸氣體積分?jǐn)?shù)最大值出現(xiàn)在臨近加熱壁面的位置。沸騰現(xiàn)象對管壁傳熱系數(shù)有明顯的影響。在壁面處的蒸氣體積分?jǐn)?shù)介于0.1～0.3之間時，壁面的傳熱系數(shù)達(dá)到最大；如果壁面處蒸氣體積分?jǐn)?shù)繼續(xù)增大，壁面?zhèn)鳠嵯禂?shù)會隨之減??；當(dāng)壁面處蒸氣體積分?jǐn)?shù)超過0.6時，壁面?zhèn)鳠嵯禂?shù)有可能會低于單相對流傳熱的傳熱系數(shù)。

符號說明：

Aq——?dú)馀菡紦?jù)的壁面面積，m2；

Alv——相間傳熱面積，m2；

a,b——系數(shù)；

cp——液相比容，W/(m·K)；

Cμb——常系數(shù)；

db——?dú)馀葜睆?，m；

dbw——?dú)馀菝撾x直徑，m；

f——?dú)馀菝撾x頻率，Hz；

fα——調(diào)節(jié)函數(shù)；

Flv——兩相間相互作用力，包括曳力、升力、湍流耗散力和壁面潤滑力，N；

h——焓，J；

hfv——蒸發(fā)潛熱，J/kg；

hl——液相對流傳熱系數(shù)，W/(m2·K)；

hlv——相間傳熱系數(shù)，W/(m2·K)；

hv——液相傳熱系數(shù)，W/(m2·K)；

hw——壁面?zhèn)鳠嵯禂?shù)，W/(m2·K)；

k——導(dǎo)熱系數(shù)，W/(m·K)；

km——混合湍流強(qiáng)度，m2/s2；

K——曳力系數(shù)；

Kbubble——?dú)馀萘饕妨粨Q系數(shù)；

Kdroplet——霧狀流曳力交換系數(shù)；

mlv——液相向氣相傳質(zhì)，代表液相的沸騰和蒸氣相的冷凝，kg；

Na——?dú)馀莩珊嗣芏龋琺-2；

p——壓力，Pa；

Qe——蒸發(fā)熱流率，W/m2；

Ql——對流傳熱熱流率，W/m2；

Qlv——相間傳熱速率，是非相變造成的相間能量傳遞，W/m3；

Qq——激冷傳熱熱流率，W/m2；

Qw——壁面總熱流率，W/m2；

Qv——蒸氣單相對流傳熱熱流率，W/m2；

r——徑向位置，m；

R0——加熱管半徑，m；

R1——玻璃管內(nèi)徑，m；

t——時間，s；

T——溫度，K；

Tsat——飽和溫度，K;

u——速度，m/s；

z——軸向，m；

α——體積分?jǐn)?shù)，%；

εm——混合湍流耗散率，m2/s2；

μ——黏性系數(shù)，Pa·s；

μeff——有效湍流黏性系數(shù)，Pa·s；

μv——?dú)庀囵ば韵禂?shù)，Pa·s；

ρ——流體密度，kg/m3；

下標(biāo):

f，w——界面，壁面；

l，v——液相，蒸氣相。

參考文獻(xiàn)

[1] BARTOLOMEI G, CHANTURIYA V. Experimental study of true void fraction when boiling subcooled water in vertical tubes[J].Thermal Engineering, 1967, 14(2):123-128.

[2] BARTEL M D. Experimental investigation of subcooled boiling[D].West Lafayette: Purdue University, 1999.

[3] LEE T H, PARK G C, LEE D J. Local flow characteristics of subcooled boiling flow of water in a vertical concentric annulus[J].International Journal of Multiphase Flow, 2002, 28(8):1351-1368.

[4] KURUL N, PODOWSKI M Z. Multidimensional effects in forced convection subcooled boiling[C]//Jerusalem： Proceedings of the Ninth International Heat Transfer Conference, 1990:19-24.

[5] ANGLART H, NYLUND O, KURUL N, et al. CFD prediction of flow and phase distribution in fuel assemblies with spacers[J].Nuclear Engineering and Design, 1997, 177(1):215-228.

[6] KONCAR B, MAVKO B. CFD simulation of subcooled flow boiling at low pressure[C]//Ljubljana：Proceedings of International Conference Nuclear Energy in Central Europe, 2001:10-13.

[7] KONCAR B, KLJENAK I, MAVKO B. Modelling of local two-phase flow parameters in upward subcooled flow boiling at low pressure[J].International Journal of Heat and Mass Transfer, 2004, 47(6):1499-1513.

[8] 李祥東， 汪榮順， 黃榮國， 等. 垂直圓管內(nèi)液氮流動沸騰的理論模型及數(shù)值模擬[J].化工學(xué)報， 2006, 57(3):491-497.(LI Xiangdong, WANG Rongshun, HUANG Rongguo, et al. Modelling and numerical simulation of boiling flow of liquid nitrogen in vertical tube[J].Journal of Chemical Industry and Engineering (China), 2006, 57(3):491-497.)

[9] LI Xiangdong, WEI Wei, WANG Rongshun, et al. Numerical and experimental investigation of heat transfer on heating surface during subcooled boiling flow of liquid nitrogen[J].International Journal of Heat and Mass Transfer, 2009, 52(5):1510-1516.

[10] 竇從從， 毛羽， 王娟， 等. 高壓高過冷度下過冷流動沸騰數(shù)值模擬[J].化工學(xué)報， 2010, 61(3):580-586.(DOU Congcong, MAO Yu, WANG Juan, et al. Numerical simulation of subcooled boiling flow under high pressure and high subcooling condition[J].Journal of Chemical Industry and Engineering (China), 2010, 61(3):580-586.)

[11] 竇從從， 毛羽， 王娟， 等. 高過冷度下直管內(nèi)氣液兩相傳熱特性[J].化工學(xué)報， 2010, 61(9):2314-2319.(DOU Congcong, MAO Yu, WANG Juan, et al. Gas-liquid two phase flow heat transfer at high subcooling in vertical tube[J].Journal of Chemical Industry and Engineering (China), 2010, 61(9):2314-2319.)

[12] MIMOUNI S, BOUCKER M, LAVIEVILLE J, et al. Modelling and computation of cavitation and boiling bubbly flows with the NEPTUNE_CFD code[J].Nuclear Engineering and Design, 2008, 238(3):680-692.

[13] UNAL H C. Maximum bubble diameter, maximum bubble growth time and bubble growth rate[J].International Journal of Heat Mass Transfer, 1976, 19(6):643-649.

[14] SATO Y, SADATOMI M, SEKOGUCHI K. Momentum and heat transfer in two-phase bubble flow—I Theory[J].International Journal of Multiphase Flow, 1981, 7(2):167-177.

[15] ISHII M, ZUBER N. Drag coefficient and relative velocity in bubbly, droplet or particulate flows[J].AIChE Journal, 1979, 25(5):843-855.

Simulation of Boiling Flow in Vertical Tube Based on Non-Equilibrium Thermal Model

FENG Liuhai1,2, NIE Yongguang1,3, WANG Jiangyun1, LEI Xiaodong4, WANG Juan1, MAO Yu1

(1.StateKeyLaboratoryofHeavyOilProcessing,ChinaUniversityofPetroleum,Beijing102249,China;2.NationalInstituteofClean-and-LowCarbonEnergy,Beijing102209,China;3.XinAoScienceandTechnologyDevelopmentCo.Ltd.,Langfang065001,China;4.ChinaMachineryKangyuanCerealsandOilEquipment(Beijing)Co.,Ltd.,Beijing100083,China)

Abstract：To simulate the boiling process, the wall heat flux was divided based on Lavievile’s thermal non-equilibrium heat flux partition algorithm. The mixture k-ε turbulent model was improved by combination of Sato model and a modifier function, by which a new drag model was presented to simulate the full boiling flow in tube. The simulated result of boiling flow in vertical circular tube agreed with experiment data well, indicating that the model was valid in simulation of boiling flow in whole tube. It is found from the simulation that the wall heat transfer coefficient was evidently enhanced in boiling, however, heat transfer deterioration appeared when the volume fraction of vapor near the wall reached about 0.6.

Key words：boiling flow; non-equilibrium thermal dynamics; drag force model; heat transfer coefficient; numerical simulation

收稿日期：2015-04-20

基金項目：國家重點(diǎn)基礎(chǔ)研究發(fā)展計劃“973”項目(2010CB226902)和中國石油大學(xué)(北京)科研基金項目(2462015YQ0303)資助

文章編號：1001-8719(2016)03-0591-06

中圖分類號：TK124

文獻(xiàn)標(biāo)識碼：A

doi:10.3969/j.issn.1001-8719.2016.03.021

第一作者： 馮留海，男，博士研究生，從事多相流動的數(shù)值模擬與實(shí)驗(yàn)方面的研究

通訊聯(lián)系人： 王娟，女，副教授，博士，從事多相流動與燃燒過程的數(shù)值模擬與實(shí)驗(yàn)研究；Tel：010-89733293；E-mail：wangjuan@cup.edu.cn；毛羽，男，教授，博士，從事多相流動及燃燒、氣固分離及液體霧化技術(shù)、化工過程裝備優(yōu)化等方面的研究；Tel：010-89733293；E-mail：maoyu@cup.edu.cn