余 健，王宏濤，廖永浩，龔光彩

（1.湖南大學(xué)土木工程學(xué)院，湖南長沙 410082； 2.91024部隊(duì)，廣東臺山 529266；3.華電電力科學(xué)研究院，浙江杭州 310030）

氣－液兩相在多孔介質(zhì)內(nèi)同向向上流動(dòng)的CFD研究＊

余 健1?，王宏濤1，2，廖永浩3，龔光彩1

（1.湖南大學(xué)土木工程學(xué)院，湖南長沙 410082； 2.91024部隊(duì)，廣東臺山 529266；3.華電電力科學(xué)研究院，浙江杭州 310030）

基于計(jì)算流體力學(xué)中歐拉模型，研究了氣－液兩相在球形多孔介質(zhì)填充床內(nèi)同向向上流動(dòng)情況.分別采用Tung ＆Siefken模型、Schulenberg模型以及Attou模型通過用戶自定義功能的形式添加氣固、液固、氣液間曳力程序?qū)μ畛浯擦鲃?dòng)情況進(jìn)行CFD三維數(shù)值模擬.對其模擬結(jié)果與實(shí)測數(shù)據(jù)進(jìn)行對比發(fā)現(xiàn)：Tung ＆Siefken的模型更適合模擬在常溫常壓下直徑在2～8 mm的球形多孔介質(zhì)內(nèi)氣－液兩相流流動(dòng)，此模型CFD數(shù)值模擬的壓降差、持液量結(jié)果均與實(shí)測數(shù)據(jù)吻合很好.

歐拉模型；氣－液兩相流；多孔材料；數(shù)值模擬

在現(xiàn)代許多工程技術(shù)領(lǐng)域中，比如在核電工業(yè)的核反應(yīng)堆、石油化工的固化床和流化床、環(huán)境工程的脫硫塔、曝氣生物濾池和滴濾池中，都存在大量關(guān)于多孔介質(zhì)內(nèi)氣－液兩相流動(dòng)問題［1－5］.因此，開展多孔介質(zhì)內(nèi)氣－液兩相流動(dòng)的研究，具有十分重要的意義.

目前，國內(nèi)外對于多孔介質(zhì)內(nèi)氣－液兩相流動(dòng)的研究多為理論分析和實(shí)驗(yàn)研究，且大部分研究都集中在氣－液同向向下流動(dòng)，而對氣－液同向向上流動(dòng)的研究較少，應(yīng)用計(jì)算流體力學(xué)（CFD）三維數(shù)值模擬氣－液同向向上流動(dòng)還尚未見報(bào)導(dǎo).Chu［1］等在試驗(yàn)的基礎(chǔ)上研究玻璃球床內(nèi)空氣和水兩相向上流動(dòng)，主要研究多孔介質(zhì)粒徑、軸線床層孔隙率的變化、氣液流量對壓降和持液量影響.Schulenberg［2］等在實(shí)驗(yàn)基礎(chǔ)上，在粗碎片多孔介質(zhì)內(nèi)考慮蒸汽與液體間曳力，提出了多孔介質(zhì)內(nèi)多相流壓降和持液量模型.Tung等［3］在多孔介質(zhì)內(nèi)兩相流實(shí)驗(yàn)的基礎(chǔ)上，提出泡狀流等不同流態(tài)概念并建立了一維壓降預(yù)測和持液量模型.Molga等［4］研究了在高壓下吸收塔內(nèi)CO2氣體與液相流動(dòng)，考慮了系統(tǒng)壓力、氣相和液相表觀速度、液相粘度以及多孔介質(zhì)粒徑對壓降和持液量影響.Schmidt［5］通過總結(jié)大多數(shù)學(xué)者關(guān)于多孔介質(zhì)內(nèi)氣液同向向上流動(dòng)的壓降模型，并在大量實(shí)驗(yàn)基礎(chǔ)上提出了自己的壓降和持液量模型，但通用性不強(qiáng).

本文首次應(yīng)用CFD三維數(shù)值模擬在常溫常壓下球形多孔介質(zhì)內(nèi)空氣和水的兩相流流動(dòng)，分別采用Tung ＆Siefken模型，Schulenberg模型以及Attou模型定性分析氣液流量和多孔介質(zhì)粒徑對壓降的影響，利用CFD對多孔介質(zhì)內(nèi)氣－液兩相同向向上流動(dòng)進(jìn)行三維數(shù)值模擬以確定最佳模型.

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

1.1 模型假設(shè)

在常溫常壓下氣－液兩相同向向上流過多孔介質(zhì)，氣體為空氣，液相為去離子水，做如下假設(shè)：

1）氣相與液相間無相變發(fā)生；

2）固定床內(nèi)的多孔介質(zhì)固定不動(dòng)；

3）球形多孔介質(zhì)被液體完全濕潤，毛細(xì)壓力可忽略；

4）球形多孔介質(zhì)被視為各向同性，不考慮軸向床層內(nèi)孔隙率的變化；

5）不考慮多孔介質(zhì)的“壁面效應(yīng)”；

6）氣體為不可壓縮流體.

1.2 體積分率連續(xù)性方程

當(dāng)氣－液兩相流流過球形多孔介質(zhì)時(shí)，每一個(gè)孔隙都由氣相與液相共同占有，即

式中：G為氣相；L為液相；ε為孔隙率.

1.3 質(zhì)量守恒方程（連續(xù)性方程）

氣相與液相流體均滿足如下質(zhì)量守恒方程：

式中：下標(biāo)K代表相，氣相K＝G，液相K＝L.

1.4 動(dòng)量守恒方程

氣相與液相流體均滿足如下動(dòng)量守恒方程：

式中：εK為K相的體積分?jǐn)?shù)；ρK為K相的密度；uK為K相的真實(shí)物理速度；P為壓力；fKR為K相與其他相動(dòng)量交換系數(shù).

1.5 NS方程封閉模型

多孔介質(zhì)內(nèi)單相流一般通過Ergun［6］公式考慮流體與多孔介質(zhì)顆粒之間的粘性和慣性作用來實(shí)現(xiàn)Navie-Stokes方程的封閉.有些學(xué)者［1－2］認(rèn)為：多孔介質(zhì)內(nèi)兩相流動(dòng)，由于一相流體的流動(dòng)影響另外一相的滲透率和通過率，建議直接對Ergun公式進(jìn)行修正以應(yīng)用于多孔介質(zhì)內(nèi)兩相流動(dòng)的研究，引入相對滲透率和相對通過率2個(gè)變量，即相對滲透率模型.不同學(xué)者提出不同的相對滲透率模型，但其模型應(yīng)用于數(shù)值模擬的可行性還需進(jìn)一步研究.

1.5.1 Tung ＆Siefken相對滲透率模型

結(jié)合單相Kozeny-Carman方程，多相流中氣固動(dòng)量交換系數(shù)的通用表達(dá)形式為：

式中：κG為氣相相對滲透率；ηG為氣相相對通過率；絕對滲透率κ和絕對通過率η采用Ergun公式［7］求得：

式中：E1和E2為Ergun系數(shù)，本文中分別取150和1.75.

Tung根據(jù)泡狀流，泡狀到彈狀過渡流，彈狀流，彈狀到柱狀過渡流以及柱狀流提出對應(yīng)流態(tài)的相對滲透率和相對通過率［3］.

式中：αL為持液量；κL為液相相對滲透率；ηL為液相相對通過率.

多孔介質(zhì)內(nèi)氣相與液相由于相對運(yùn)動(dòng)從而產(chǎn)生相互作用曳力，其值的確定在國際上一直是難點(diǎn).Schulenberg［2］，Tung［3］，Schmidt［5］，Attou［8］等提出了不同的氣液曳力模型.Tung認(rèn)為當(dāng)流態(tài)為泡狀流和彈狀流時(shí)，氣相以單一分散氣泡的形式在多孔介質(zhì)內(nèi)運(yùn)動(dòng)；在較低流速（Stokes流）時(shí)，粘性力在氣泡與液相間曳力中占主導(dǎo)地位；在較高流速時(shí)，曳力主要來自慣性力.綜合考慮粘性力與慣性力，提出下面的曳力系數(shù)表達(dá)式［3］：

Tung在試驗(yàn)的基礎(chǔ)上，根據(jù)泡狀流和彈狀流提出對應(yīng)流態(tài)的粘性系數(shù)C和慣性系數(shù)C［3］.

式中：γ為氣泡與顆粒直徑比；σ為表面張力；f為形狀影響因子；dP為多孔介質(zhì)當(dāng)量粒徑.

關(guān)于彈狀到柱狀過渡流和柱狀流，Tung認(rèn)為彈狀氣泡合并在一起形成氣相連續(xù)相在通道流動(dòng)，而液相以液滴的形式流動(dòng)［3］，其氣液曳力系數(shù)表達(dá)式分別為式（28）和式（29）.

式中：κG和ηG采用式（7）計(jì)算；W采用公式（8）計(jì)算.

式中：κG和ηG采用公式（9）計(jì)算.

1.5.2 Schulenberg模型

Schulenberg基于兩相流的動(dòng)量方程研究多孔介質(zhì)內(nèi)氣液流動(dòng)，綜合考慮氣液間的虛擬質(zhì)量力，毛細(xì)力，粘性力與慣性力對曳力影響，提出氣固、液固和氣液曳力模型［5］.其氣固、液固和氣液間動(dòng)量交換系數(shù)分別為式（30），式（33）和式（35）.

1.5.3 Attou雙流體界面模型

Attou忽略層流粘性應(yīng)力張量，用面積平均法建立穩(wěn)態(tài)滴濾池一維氣液兩相數(shù)學(xué)模型，提出氣固、液固和氣液曳力模型［8］.氣固、液固和氣液間動(dòng)量交換系數(shù)為：

2 模擬對象與邊界條件

2.1 模擬對象

試驗(yàn)裝置如圖1（a）所示.采用內(nèi)徑為50 mm，長度為1 m的有機(jī)玻璃管，豎直固定放置，管中分別填充直徑為2，5，8 mm的玻璃球，測壓環(huán)之間的距離為700 mm，氣體（空氣）與液體（去離子水）經(jīng)氣液混合器從下面流入球形多孔介質(zhì)球床，從裝置上端流出.氣體與液體流量通過流量計(jì)測量，測壓環(huán)的壓強(qiáng)通過差壓變送器測得，液體質(zhì)量流量為36～324 kg／h，氣體質(zhì)量流量為0.05～13 kg／h［9］.

本文的模擬對裝置進(jìn)行合理簡化，反應(yīng)器形狀、尺寸如圖1（b）所示.模擬計(jì)算區(qū)域?yàn)椋簻y壓環(huán)之間距離700 mm的試驗(yàn)段，其流體物質(zhì)屬性與實(shí)測值相同.利用ICEM對計(jì)算區(qū)域進(jìn)行三維網(wǎng)格劃分，網(wǎng)格類型為0型網(wǎng)格，網(wǎng)格單元總數(shù)為30 720個(gè).

圖1 試驗(yàn)裝置與CFD物理模型Fig.1 Diagram of experiment and CFD physical model

2.2 模型設(shè)置

選用Fluent 6.3軟件中多相流模型，即氣液固三相歐拉模型，氣體為連續(xù)相phase1，液體為phase2，球形多孔介質(zhì)為靜止相phase3，氣相、液相和多孔介質(zhì)三相間的動(dòng)量交換系數(shù)通過編寫用戶自定義功能（UDF）添加到程序中.采用壓力基非穩(wěn)態(tài)求解器，時(shí)間步長為0.01 s，方程離散格式為一階迎風(fēng)格式，壓力和速度的耦合采用simple算法，質(zhì)量殘差設(shè)定為5×10－4，其他殘差設(shè)置為2×10－6.

2.3 初始條件和邊界條件

反應(yīng)器入口為速度入口，按照實(shí)驗(yàn)進(jìn)口處測定氣相與液相的質(zhì)量流量分別轉(zhuǎn)化為對應(yīng)的速度，設(shè)定液相入口的體積分?jǐn)?shù)，出口均為自由出流邊界條件，壁面采用無滑移邊界條件.

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

3.1 阻力壓降

影響氣液兩相流過多孔介質(zhì)壓降的因素有很多，如多孔介質(zhì)顆粒的粒徑、氣液兩相的流量、流體的屬性和孔隙率的大小等眾多因素.本文基于Tung ＆Siefken模型、Schulenberg模型以及Attou模型，利用CFD數(shù)值模擬方法，研究了氣液流量、孔隙率和球形多孔介質(zhì)粒徑對阻力壓降的影響，結(jié)果如圖2，圖3和圖4所示.

圖2 進(jìn)液量為0.045 kg／s的2 mm球床各模型壓降與試驗(yàn)數(shù)據(jù)對比Fig.2 Comparison between pressure-drop models and experimental results（dp＝2mm，L＝0.045 kg／s）

從圖2～圖4可知，當(dāng)固定進(jìn)液流量，進(jìn)氣流量較小時(shí)，阻力壓降增加較緩慢；當(dāng)進(jìn)氣量逐漸增大時(shí)，阻力壓降增加較快，這主要是由于阻力壓降包括粘性力分量和慣性力分量，粘性力與速度成比例變化關(guān)系，慣性力與速度的平方成比例變化關(guān)系所致［10］.當(dāng)進(jìn)氣流速較小時(shí)，氣相的雷諾數(shù)較小，粘性力分量所占份額較大，從而阻力壓降增加較緩慢.隨著進(jìn)氣流速增加，其雷諾數(shù)較高，氣相的慣性力所占份額較大，導(dǎo)致氣體與多孔介質(zhì)邊壁、氣液間的作用力越來越大，并且較高的氣相雷諾數(shù)推動(dòng)同向的液體流動(dòng)［3］，從而導(dǎo)致阻力壓降增加較快.

通過模擬和試驗(yàn)均發(fā)現(xiàn)，同樣的氣體和液體流量，顆粒越小阻力越大，這主要由于顆粒越小，其孔隙率越小，氣液兩相流在多孔介質(zhì)內(nèi)的真實(shí)物理速度越大，從而導(dǎo)致氣體與多孔介質(zhì)邊壁、液體與多孔介質(zhì)邊壁以及氣液間的作用劇烈.

圖3 進(jìn)液量為0.031 kg／s的5 mm球床各模型壓降與試驗(yàn)數(shù)據(jù)對比Fig.3 Comparison between pressure-drop models and experimental results（dp＝5 mm，L＝0.031 kg／s）

圖4 進(jìn)液量為0.088 kg／s的8 mm球床各模型壓降與試驗(yàn)數(shù)據(jù)對比Fig.4 Comparison between pressure-drop models and experimental results（dp＝8 mm，L＝0.088 kg／s）

通過圖2，圖3和圖4發(fā)現(xiàn)，Schulenberg模型模擬值均遠(yuǎn)大于實(shí)測值，原因可能是：其模型主要應(yīng)用于高溫高壓下核反應(yīng)堆多孔介質(zhì)內(nèi)氣液兩相流動(dòng)，其計(jì)算的液體相對通過率遠(yuǎn)小于其他模型的相對滲透率［2］，導(dǎo)致計(jì)算的液相慣性力較大，從而使模擬值遠(yuǎn)大于常溫常壓下的實(shí)測值.

Attou模擬值均比實(shí)測值偏大，其原因可能是：Attou模型是基于環(huán)狀流流態(tài)基礎(chǔ)上建立的雙流體界面數(shù)學(xué)模型，主要應(yīng)用于滴濾池內(nèi)氣液兩相向下流動(dòng)［8］，其模型應(yīng)用于多孔介質(zhì)內(nèi)非柱狀流動(dòng)還值得商榷.通過數(shù)值模擬發(fā)現(xiàn)，基于Attou模型不適合應(yīng)用于多孔介質(zhì)內(nèi)氣液兩相同向向上流動(dòng).

在常溫常壓下濾料直徑在2～8 mm的多孔介質(zhì)床內(nèi)氣液兩相流流動(dòng)中，利用Tung ＆Siefken模型數(shù)值模擬壓降值均與實(shí)測結(jié)果吻合很好，其主要原因可能是：Tung ＆Siefken模型考慮了泡狀流，泡狀－彈狀流混合流，彈狀流，彈狀－環(huán)狀混合流以及環(huán)狀流不同的流型，并建立對應(yīng)的曳力計(jì)算公式.

圖5為利用Tung ＆Siefken模型多孔介質(zhì)濾料粒徑為8 mm，ReL為272.8，ReG為605.6數(shù)值模擬壓降沿裝置軸面變化云圖.通過模擬可知，壓降沿著流動(dòng)方向均勻降低，且與實(shí)測壓降數(shù)據(jù)較好吻合.

圖5 dp＝8 mm，ReL＝272.8，ReG＝605.6軸面壓降云圖Fig.5 The contours of pressure-drop changes along axial plane（dp＝8 mm，ReL＝272.8，ReG＝605.6）

圖6為采用Tung ＆Siefken模型液體質(zhì)量流量為36～324 kg／h，氣體質(zhì)量流量為0.05～13 kg／h，進(jìn)行CFD數(shù)值模擬值與實(shí)測值的關(guān)系圖，并采用如下統(tǒng)計(jì)公式進(jìn)行誤差分析：

圖6 Tung ＆Siefken模型CFD模擬壓降值與實(shí)測值對比Fig.6 Comparison of pressure-drop based on Tung ＆Siefken CFD simulation with experimental results

平均相對誤差ē為7.3%和相對誤差的均方根erms為10.2%，通過圖6與誤差統(tǒng)計(jì)分析均得出：Tung ＆Siefken的模型能較準(zhǔn)確數(shù)值模擬粒徑為2～8 mm的多孔介質(zhì)球床內(nèi)氣液兩相流的壓降.

3.2 持液量

持液量為多孔介質(zhì)空隙內(nèi)液體相的體積分?jǐn)?shù)，不同的進(jìn)氣量和進(jìn)液量以及孔隙率直接影響持液量值.圖7顯示了在液體流量為0.065 kg／s時(shí)，基于Tung ＆Siefken模型的CFD模擬持液量值隨進(jìn)氣量與顆粒直徑的變化情況.

圖7 進(jìn)液量為0.065 kg／s的CFD數(shù)值模擬持液量與試驗(yàn)數(shù)據(jù)對比Fig.7 Comparison of liquid holdups using CFD simulation with experimental results（L＝0.065 kg／s，dP＝2 mm，dP＝5 mm）

由圖7可知，CFD數(shù)值計(jì)算持液量基本與實(shí)測數(shù)據(jù)一致，持液量隨進(jìn)氣量的增加而逐漸降低.當(dāng)粒徑為2 mm時(shí)，液體持液量的CFD模擬值均大于實(shí)測值；當(dāng)粒徑為5 mm時(shí)，液體持液量的大多數(shù)CFD模擬值小于實(shí)測值.出現(xiàn)誤差的原因是由于Tung＆Siefken模型是在關(guān)于阻力壓降實(shí)驗(yàn)的基礎(chǔ)上總結(jié)的半經(jīng)驗(yàn)?zāi)Ｐ?，利用此模型?jì)算出的持液量有一定誤差，故需要通過大量實(shí)驗(yàn)對模型進(jìn)行進(jìn)一步完善.

4 結(jié) 論

在常溫常壓ReL＝33～300，ReG＝10～650內(nèi)，本文基于Fluent計(jì)算流體力學(xué)軟件中歐拉模型，分別采用Tung ＆Siefken模型、Schulenberg模型以及Attou模型，通過用戶自定義功能編寫氣固、液固、氣液間的曳力系數(shù)程序到軟件中的形式，首次對球形多孔介質(zhì)填充床內(nèi)氣液兩相向上流動(dòng)進(jìn)行CFD三維數(shù)值模擬.通過模擬發(fā)現(xiàn)，Tung ＆Siefken模型更適合應(yīng)用于粒徑為2～8 mm的多孔介質(zhì)球體內(nèi)氣液兩相流的數(shù)值模擬，模擬壓降值與實(shí)測值基本一致，平均相對誤差為7.3%.CFD數(shù)值模擬的持液量與實(shí)測值變化基本一致，但是存在一定的誤差.

［1］ CHU W，DHIR V K，MARSHALL J.Study of pressure drop，void fraction and relative permeabilities of two phase flow through porous media［J］.Heat Transfer-Seattle，AIChE Symposium series，1983，79：224－235.

［2］ SCHULENBERG T，MULLER U.A refined model for the coolability of core debris with flow entry from the bottom［C］／／Proceedings of the Sixth Information Exchange Meeting on Debris Coolability.Los Angeles，1984：1090－1097.

［3］ TUNG V X，DHIR V K.A hydrodynamic model for twophase flow through porous media［J］.International Journal of Multiphase Flow，1988，14：47－65.

［4］ MOLGA E J，WESTERTERP K R.Gas-liquid interfacial area and holdup in a cocurrent upflow packed bed bubble column reactor at elevated pressures［J］.Industrial ＆Engineering Chemistry Research，1997，36：622－631.

［5］ SCHMIDT W.Interfacial drag of two-phase flow in porous media［J］.International Journal of Multiphase Flow，2007，33：638－657.

［6］ ERGUN S.Fluid flow through packed columns［J］.Chemical Engineering Progress，1952，48：89－105.

［7］ SIEFKEN L J，CORYELL E W.SCDAP／RELAP5 modeling of heat transfer and flow losses in lower head porous debris［R］.Idaho Falls：Idaho National Laboratory，1999.

［8］ ATTOU A，F(xiàn)ERSCHNEIDERR G.A two-fluid hydrodynamic model for the transition between trickle and pulse flow in a cocurrent gas-liquid packed bed reactor［J］.Chemical Engineering Science，2000，55：491－511.

［9］ 李振鵬.球床多孔介質(zhì)通道單相流體流動(dòng)特性研究［D］.哈爾濱：哈爾濱工程大學(xué)核科學(xué)與技術(shù)學(xué)院，2009.

LI Zhen-peng.Study on single phase flow in pebble-bed porous media channel［D］.Harbin：College of Nuclear Science and Technology，Harbin Engineering University，2009.（In Chinese）

［10］廖永浩.微球床多孔介質(zhì)通道氣液兩流特性研究［D］.哈爾濱：哈爾濱工程大學(xué)核科學(xué)與技術(shù)學(xué)院，2009.

LIAO Yong-hao.Study on gas-1iquid two phase flow in microspheres bed porous media channel［D］.Harbin：College of Nuclear Science and Technology，Harbin Engineering University，2009.（In Chinese）

CFD Simulation of Gas-liquid Two Phase Cocurrent Upflow through Porous Media

YU Jian1?，WANG Hong-tao1，2，LIAO Yong-hao3，GONG Guang-cai1

（1.College of Civil Engineering，Hunan Univ，Changsha，Hunan 410082，China；2.Army 91024，PLA，Taishan，Guangdong 529266，China；
3.Huadian Electric Power Research Institute，Hangzhou，Zhejiang 310030，China）

Based on computational fluid dynamics Euler multiphase model，gas-liquid two-phase co-current upflow through packed bed was investigated.Tung ＆Siefken model，Schulenberg model and Attou k-fluid model were applied，respectively.The gas-particle，liquid-particle and gas-liquid drag force programs were implemented to CFD three-dimensional numerical simulation with user defined functions.A comparison of the experiment data has shown that Tung ＆Siefken＇s model was more suitable for the numerical simulation of gas-liquid two-phase upflow through porous media packed with 2～8 mm sphere under normal temperature and pressure.The CFD numerical simulated values of pressure drop and liquid hold-up by Tung ＆Siefken＇s model were consistent with the experiment results.

Euler model；gas-liquid flow；porous materials；numerical simulation

X11

A

1674-2974（2012）08-0067-06＊

2011-11-21

國家“十一五”水體污染控制與治理科技重大專項(xiàng)（2009ZX07423－004）

余 ?。?964—），男，湖南益陽人，湖南大學(xué)教授，博士生導(dǎo)師

?通訊聯(lián)系人，E-mail：jianyu＠hnu.edu.cn