收稿日期： 2022-03-23； 修回日期： 2022-05-25； 網(wǎng)絡(luò)出版時間： 2022-09-22

網(wǎng)絡(luò)出版地址： https：//link.cnki.net/urlid/32.1814.TH.20220919.1832.009

基金項目： 國家自然科學(xué)基金資助項目（52106216）

第一作者簡介： 李玥嬛（1998—），女，內(nèi)蒙古通遼人，碩士研究生（1808320994@qq.com），主要從事氣固兩相流數(shù)值模擬研究.

通信作者簡介： 朱曉麗（1991—），女，山東濰坊人，講師，博士（zhuxiaoli@upc.edu.cn），主要從事氣固兩相流數(shù)值模擬研究.

摘要： 為明確加壓鼓泡流化床內(nèi)氣固兩相流動基本規(guī)律以及壓力的影響機制，為多種工業(yè)過程領(lǐng)域加壓鼓泡流化床反應(yīng)器的設(shè)計、運行和放大提供有益參考.文中基于CFD-DEM方法，以開源軟件Mfix為計算框架，開展了加壓鼓泡流化床氣固兩相流動數(shù)值模擬研究.首先依托薄矩形加壓鼓泡流化床試驗臺數(shù)據(jù)進行模型驗證，在此基礎(chǔ)上，系統(tǒng)地研究了操作壓力對Geldart B類顆粒臨界流化速度的影響，獲得了不同操作壓力下床層內(nèi)顆粒運動與分布規(guī)律以及氣泡尺寸和頻率分布.數(shù)值模擬結(jié)果表明：在0.1～0.7 MPa，隨著操作壓力的升高，Geldart B類顆粒臨界流化速度減小，并且在低操作壓力下臨界流化速度隨操作壓力變化更加明顯；同時，隨著操作壓力的升高，床層內(nèi)氣泡尺寸減小，床層均勻性有所改善.

關(guān)鍵詞： 加壓鼓泡床；氣固流動特性；臨界流化速度；氣泡特性

中圖分類號： S277.9" 文獻標(biāo)志碼： A" 文章編號： 1674-8530（2024）06-0570-06

DOI：10.3969/j.issn.1674-8530.22.0069

李玥嬛，朱曉麗，王振波，等.基于CFD-DEM方法的加壓鼓泡床氣固流動特性數(shù)值模擬［J］.排灌機械工程學(xué)報，2024，42（6）：570-575.

LI Yuehuan，ZHU Xiaoli，WANG Zhenbo， et al. Numerical simulation of gas-solid flow characteristics in pressurized bubbling flui-dized bed based on CFD-DEM method［J］.Journal of drainage and irrigation machinery engineering（JDIME），2024，42（6）：570-575.（in Chinese）

Numerical simulation of gas-solid flow characteristics in pressurized

bubbling fluidized bed based on CFD-DEM method

LI Yuehuan， ZHU Xiaoli*， WANG Zhenbo， LIU Yibo

（College of New Energy， China University of Petroleum （East China）， Qingdao， Shandong 266580， China）

Abstract： Numerical simulation study on a pressurized bubbling fluidized bed was carried out based on a coupled CFD-DEM approach using an open-source software （Mfix） in order to clarify the basic laws of gas-solid two-phase flow in a pressurized bubbling fluidized bed and the influencing mechanism of pressure， which provides useful reference for the design， operation and scale-up of pressurized bubbling fluidized bed reactors in various industrial processes. Firstly， the numerical model was verified by the experimental data from the thin rectangular pressurized bubbling fluidized bed. Based on this， the influence of operating pressure on the critical fluidization velocity of Geldart B-type particles， the particle movement and distribution laws， as well as the bubble size and frequency distribution at different operating pressures were systematically studied. The numerical simulation results show that the critical fluidization velocity of Geldart B-type particles decreases with the increase of operating pressure in the range of 0.1 MPa to 0.7 MPa， and the critical fluidization velocity changes more obviously with the increasing pressure under low operating pressure conditions. In addition， as pressure increases， the bubble size decreases， leading to improvement in bed uniformity.

Key words： pressurized bubbling fluidized bed;gas-solid flow characteristics;critical fluidization velocity;bubble characteristics

鼓泡流化床反應(yīng)器（以下簡稱鼓泡床）廣泛應(yīng)用于煤炭燃燒與氣化、冶金工業(yè)、制藥裝備等領(lǐng)域.高壓操作可有效促進氣固混合、強化傳質(zhì)傳熱以及減小設(shè)備尺寸［1］，是工業(yè)鼓泡流化床反應(yīng)器的重要發(fā)展趨勢.多年來，研究人員對于常壓鼓泡床氣固流動特性已經(jīng)有了較為充分的認(rèn)識和理解，然而在加壓條件下，床層內(nèi)氣泡的尺寸、頻率以及上升速度等基本流動特性都會發(fā)生改變，進而影響鼓泡流化床反應(yīng)器的整體性能，因此有必要對加壓鼓泡床氣固流動特性開展研究.

針對加壓鼓泡流化床，前人已開展大量試驗研究工作.然而，目前加壓條件下氣泡的基本流動特征規(guī)律在文獻中并未形成統(tǒng)一的結(jié)論.部分研究人員如SONG等［2］、OLOWSON等［3-4］以及GODLIEB等［5-6］發(fā)現(xiàn)，隨著壓力的升高，氣泡尺寸以及氣泡上升速度均有所減小.HOFFMANN等［7］認(rèn)為，氣泡在床層內(nèi)所處的位置以及操作壓力的范圍等都會影響氣泡尺寸.

對于氣固流態(tài)化體系，試驗測量雖然能夠提供更為準(zhǔn)確的流場參數(shù)，但試驗測量很難反映加壓鼓泡床內(nèi)復(fù)雜的多尺度氣固流動信息，特別是顆粒尺度上的細節(jié)信息.計算流體力學(xué)（CFD）的發(fā)展能夠彌補現(xiàn)有試驗條件和測試手段的不足.常用的描述氣固兩相流動的數(shù)值方法包括基于歐拉-歐拉框架的TFM方法［8］以及基于歐拉-拉格朗日框架的MP-PIC［9-10］和CFD-DEM方法［11-13］.其中，CFD-DEM方法可追蹤真實的單顆粒的運動與受力，具有很高的流場信息分辨率.MANSOURPOUR等［14］和LIU等［15］采用CFD-DEM方法對加壓鼓泡流化床氣固兩相流動開展研究，得到氣泡尺寸和氣泡速度等隨壓力的變化關(guān)系.宋加龍等［16］發(fā)現(xiàn)了操作壓力和流化數(shù)對氣泡特性包括氣泡尺寸、上升速度與氣泡體積分?jǐn)?shù)以及顆?；旌系确矫娴挠绊?

綜上，目前加壓鼓泡流化床內(nèi)氣固流動特性的基本規(guī)律以及壓力的影響機制尚不明確.所以，文中采用CFD-DEM方法分析壓力對臨界流化速度、顆粒速度分布以及氣泡尺寸和頻率分布的影響，以期為加壓鼓泡流化床反應(yīng)器的運行和設(shè)計提供一定的理論參考.

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

1.1" 氣相控制方程

氣相質(zhì)量守恒和動量守恒方程為

t（ρgεg）+SymbolQC@ ·（ρgεgug）=0，（1）

t（ρgεgug）+SymbolQC@ ·（ρgεgugug）=

－εgSymbolQC@ pg+

SymbolQC@ ·εgτg+εgρgg－F，（2）

式中：εg為空隙率，即網(wǎng)格內(nèi)的氣相體積分?jǐn)?shù);ρg為氣體密度；ug為速度矢量;pg為氣體壓力;g為重力加速度;τg為氣體黏性應(yīng)力張量；F為顆粒相與氣相之間的動量交換源項.

1.2" 固相控制方程

顆粒相的運動遵循牛頓第二定律，作用于顆粒上的力主要包括壓力梯度力、重力、氣固相間曳力以及顆粒間碰撞作用力等.對于單顆粒a，其平動和轉(zhuǎn)動方程可以描述為

madvadt=－Va

SymbolQC@ pg+mag+Fd+Fc，（3）

Iadωadt=∑N1b=1（Rbnab×Ftab），（4）

式中：ma，Va分別為顆粒a的質(zhì)量和體積；va，ωa為顆粒平動速度以及角速度；Ia為轉(zhuǎn)動慣量；Fc為顆粒間碰撞作用力；Fd為顆粒所受曳力，F(xiàn)d=Vaβεp（ug－up），εp為當(dāng)前計算網(wǎng)格內(nèi)的顆粒相體積分?jǐn)?shù)，ug為當(dāng)前網(wǎng)格內(nèi)流體速度矢量，up為顆粒速度矢量，β為相間曳力系數(shù)；N1為與顆粒a接觸作用的顆?？倲?shù)；Rb為顆粒之間的距離；nab為從顆粒a指向顆粒b的法向單位矢量；Ftab為顆粒切向碰撞力.

2" 試驗及模擬參數(shù)設(shè)置

自行搭建尺寸為100 mm×34 mm×500 mm的加壓鼓泡流化床試驗臺，由鼓泡床本體、供風(fēng)系統(tǒng)以及測量系統(tǒng)3部分構(gòu)成.對Geldart B類顆粒開展試驗和數(shù)值模擬研究，通過壓差傳感器測量獲得床層壓降-氣速曲線，進而得到臨界流化速度，從而對數(shù)值模型的準(zhǔn)確性進行驗證.試驗中所采用的玻璃珠顆粒具有寬篩分粒徑分布（410～1 040 μm，近似正態(tài)分布），其平均粒徑為721 μm，密度為2 550 kg/m3，為簡化計算，將該平均值作為均一粒徑開展數(shù)值模擬.顆粒初始堆積高度為0.08 m，操作壓力p分別為0.1，0.4，0.7 MPa，對應(yīng)氣相密度分別為1.21，4.84和8.47 kg/m3，氣相黏度為1.8×10-5 kg/（m

瘙 簚 s），流化數(shù)（表觀氣速U/臨界流化速度Umf）為3.設(shè)置摩擦系數(shù)為0.5，法向及切向恢復(fù)系數(shù)為0.9.試驗臺及邊界條件設(shè)置如圖1所示，圖中H0為初始床層高度.

為減小計算量和降低時間成本，對模型進行適當(dāng)簡化，在Z方向僅設(shè)置3層網(wǎng)格，并設(shè)置為周期性邊界［17］，減小模型尺寸及計算域內(nèi)的網(wǎng)格和顆粒數(shù)量，從而提高計算效率.每個工況計算5 s，每秒保存80組數(shù)據(jù)，為避免啟動過程對數(shù)據(jù)結(jié)果產(chǎn)生影響，取后3 s準(zhǔn)穩(wěn)態(tài)階段數(shù)據(jù)進行分析.

3" 氣泡圖像的處理方法

為提取床層內(nèi)氣泡形狀、尺寸與空間信息，借助MATLAB軟件，對數(shù)值模擬圖像進行處理.首先去除圖像背景，并確保被裁剪后的圖像長寬值固定，將圖像依次進行灰度變換和二值轉(zhuǎn)換，將灰度圖像中像素值大于閾值（采用氣相體積分?jǐn)?shù)0.8作為閾值，符號表示為Th）的設(shè)置為1，認(rèn)為是氣泡相，小于Th的設(shè)置為0，認(rèn)為是乳化相.為提高數(shù)據(jù)的準(zhǔn)確性，將二值化之后的圖像進行均值濾波、區(qū)域聯(lián)通等形態(tài)學(xué)處理.如圖2所示，數(shù)值模擬所得到的瞬態(tài)顆粒分布圖像，經(jīng)過圖像處理后可提取到床層內(nèi)的5個氣泡.

進一步對二值圖像進行統(tǒng)計，計算某一氣泡所占的像素點數(shù)量乘以每個像素點對應(yīng)的實際面積作為氣泡面積，進而求解氣泡當(dāng)量直徑.具體計算式為

Ab=∑b（Xj，Yj）=∮b（x，y）dxdy，（5）

Db=4Abπ，（6）

式中：Ab為氣泡投影面積；b（Xj，Yj）為第j個像素點的面積；b（x，y）為氣泡邊界曲線方程；Db為氣泡當(dāng)量直徑.

為探究氣泡在床層內(nèi)的時空分布規(guī)律，每種計算工況下，提取計算時間3 s內(nèi)共240幅采集圖像中所有氣泡進行數(shù)據(jù)統(tǒng)計.根據(jù)氣泡的質(zhì)心坐標(biāo)統(tǒng)計其軸向高度的氣泡數(shù)量，將某一軸向高度的氣泡數(shù)量除以氣泡總數(shù)量作為該軸向高度的氣泡頻率.計算方法為

f=Ni，jN2，（7）

式中： f為氣泡頻率；Ni，j為空間位置為（i，j）處的氣泡個數(shù)；N2為床層內(nèi)的氣泡總數(shù).

4" 結(jié)果及分析

4.1" 臨界流化速度

一般而言，臨界流化速度Umf根據(jù)床層壓降與氣速的關(guān)系確定.在固定床階段，壓降隨氣速增加近似呈線性增加，當(dāng)床層進入流化之后其壓降基本保持不變，固定床壓降曲線與流化床壓降曲線的交點即對應(yīng)Umf.在流態(tài)化階段，固體顆粒所受到的曳力等于顆粒重力，顆粒處于懸浮狀態(tài)，如果不考慮流體和顆粒與壁面之間的摩擦力，則床層壓降全部轉(zhuǎn)化為流體對顆粒的曳力，即在理想狀態(tài)下，床層壓降等于單位面積內(nèi)顆粒的重量，故采用單位面積的顆粒重量對床層壓降進行歸一化.圖3為不同壓力下，歸一化床層壓降Δpn與表觀氣速U的關(guān)系曲線.

從圖3可以看出，在不同操作壓力下，流態(tài)化階段的歸一化壓降值在1.0～1.1，略大于理想狀態(tài)時的歸一化壓降值（Δpn=1.0），多余的壓力可能與壁面效應(yīng)和顆粒之間的摩擦等有關(guān).當(dāng)壓力從0.1 MPa增加到0.4 MPa，數(shù)值模擬中Umf從0.480 m/s降至0.280 m/s，試驗測量中Umf從0.405 m/s降至0.270 m/s；當(dāng)壓力從0.4 MPa增加到0.7 MPa時，數(shù)值模擬中Umf從0.280 m/s降至0.225 m/s，試驗測量中Umf從0.270 m/s降至0.210 m/s.可以看出，在較低的操作壓力下（0.1～0.4 MPa），臨界流化速度隨壓力的變化更為明顯.

圖4為顆粒臨界流化速度模擬值和試驗值對比.

由圖4可見，不同壓力下數(shù)值模擬與試驗結(jié)果誤差均在20%以內(nèi)，誤差在合理范圍內(nèi).模擬值與試驗值差異的原因可能與曳力模型的選取有關(guān)，文中采用的Syamlal-O′Brien曳力模型基于均勻流態(tài)化假設(shè)，未能考慮到介尺度氣泡和顆粒團聚效應(yīng)對曳力的影響，因此導(dǎo)致試驗結(jié)果與數(shù)值模擬有一定偏差.

操作壓力升高，臨界流化速度減小這一現(xiàn)象可通過Ergun方程［18］進行解釋，即

ΔpH=150（1-εmf）2ε3mfμUmf（sdp）2+1.751－εmfε3mfρgU2mfsdp，（8）

式中：Δp為壓降；H為床層軸向高度;dp為顆粒直徑;s為顆粒球形度（球形顆粒s=1）；εmf為臨界流化速度下的床層孔隙率；μ為氣體黏度.方程左側(cè)為床層壓降，方程右側(cè)第一項為黏性力所引起的壓力耗散，第二項為慣性力所引起的壓力耗散.當(dāng)操作壓力升高時，氣體密度增大，而氣體黏度幾乎不受影響，因此為保證恒定的床層壓降，臨界流化速度則相應(yīng)減小.

4.2" 顆粒速度分布

圖5為不同操作壓力下顆?？臻g位置分布，圖中vy是顆粒的軸向速度.從圖中可以看出，當(dāng)操作壓力從0.1 MPa升高至0.7 MPa時，床層內(nèi)氣泡尺寸減小，氣泡數(shù)量增加，乳化相在床層內(nèi)分布更加均勻.

為分析壓力對床層內(nèi)顆粒空間分布和運動特性的影響，提取床層內(nèi)不同空間位置處顆粒軸向時均速度，做出不同床層軸向高度處顆粒軸向時均速度vs沿床層水平方向的分布曲線，3種操作壓力下的對比如圖6所示，圖中x/X為相對徑向位置，X為床體寬，x是選取的測量位置.從圖6中可以看出，顆粒軸向時均速度呈軸對稱分布，床層中心區(qū)域形成顆粒軸向時均速度為正的上行流，邊壁區(qū)域形成顆粒時均速度為負(fù)的下行流.從圖6a中可以看出，在常壓下，隨著床層軸向高度升高，核心區(qū)域逐漸變大，這是由于在較高的床層軸向高度處小氣泡聚集并成大氣泡，床層不均勻性增強，因此核心區(qū)域逐漸增大.隨著壓力的升高，核心區(qū)域基本不隨床層軸向高度變化，但核心區(qū)域整體向壁面方向擴展.這是由于隨著壓力升高，氣泡尺寸減小，氣泡數(shù)量增加，氣泡在徑向分布更加均勻，因此導(dǎo)致核心區(qū)域增大.

4.3" 氣泡尺寸及頻率分布

圖7為不同操作壓力下瞬態(tài)氣相體積分?jǐn)?shù)φ分布，從圖中可以看出，隨著壓力升高，氣泡數(shù)量增多，氣泡尺寸明顯減小.

圖8為氣泡當(dāng)量直徑在不同操作壓力下隨床層軸向高度的變化曲線.

如圖8所示，為進一步探究氣泡尺寸與壓力的關(guān)系，沿床層軸向由高到低選取A1—A6這6個軸向高度位置（見圖1）為監(jiān)測點，隨著軸向高度增加，氣泡當(dāng)量直徑增大.這是由于隨著軸向高度增大，床層壓力降低，氣泡在上升過程中逐漸變大并且相互聚集.隨操作壓力的升高，氣泡當(dāng)量直徑明顯減小，與圖7中展示的瞬態(tài)氣相體積分?jǐn)?shù)分布規(guī)律一致，且與前人的研究結(jié)果一致［2］.因為操作壓力的升高加強了氣泡間相互作用，使氣泡穩(wěn)定性降低，氣泡傾向于分裂成更小的單元，從而導(dǎo)致氣泡尺寸減小.

圖9為不同操作壓力下氣泡頻率軸向分布曲線.

從圖9可以看出，在不同操作壓力下氣泡頻率沿軸向均呈單峰分布，且峰值隨操作壓力增大而增大，這是由于氣泡尺寸隨操作壓力增大而減小，因此導(dǎo)致氣泡頻率的增加.隨著床層軸向高度增加，氣泡頻率呈顯著減小趨勢，表明在該床段氣泡以聚并為主.隨著床層軸向高度進一步增加，在常壓下氣泡頻率隨軸向高度增加先趨于穩(wěn)定后逐漸下降，表明氣泡沿軸向高度穩(wěn)定上升并伴隨生長聚并，直到達到床層表面發(fā)生破裂.而在高壓下，隨著軸向高度的增加，氣泡頻率呈現(xiàn)波動下降狀態(tài)，這表明在高壓下氣泡聚并破碎程度加劇，且氣泡破碎現(xiàn)象和常壓狀態(tài)相比明顯增多.

5 "結(jié)" 論

1） 隨著操作壓力升高，顆粒臨界流化速度減小，且低操作壓力下臨界流化速度的變化更加顯著.

2） 隨著操作壓力升高，核心區(qū)域擴散程度增加，床層均勻性提升.

3） 隨著操作壓力升高，氣泡穩(wěn)定性降低，氣泡當(dāng)量直徑減小，氣泡數(shù)量增加.

參考文獻（References）

［1］" 朱治平. 加壓循環(huán)流化床的實驗與模型研究［D］. 北京：中國科學(xué)院工程熱物理研究所， 2008.

［2］" SONG J， LIU D， MA J， et al. Effect of elevated pressure on bubble properties in a two-dimensional gas-solid fluidized bed［J］. Chemical engineering research amp; design， 2018， 138：21-31.

［3］" OLOWSON P A， ALMSTEDT A E. Influence of pressure and fluidization velocity on the bubble behaviour and gas flow distribution in a fluidized bed［J］. Chemical engineering science， 1990， 45（7）：1733-1741.

［4］" OLOWSON P A， ALMSTEDT A E. Hydrodynamics of a bubbling fluidized bed： influence of pressure and flui-dization velocity in terms of drag force［J］. Chemical engineering science， 1992， 47（2）：357-366.

［5］" GODLIEB W， DEEN N G， KUIPERS J A M. Bubble behaviour in fluidised beds at elevated pressures［J］. Macromolecular materials and engineering， 2011， 296（3/4）：270-277.

［6］" GODLIEB W， GORTER S， DEEN N G， et al. Experimental study of large scale fluidized beds at elevated pressure［J］. Industrial amp; engineering chemistry research， 2012， 51（4）：1962-1969.

［7］" HOFFMANN A C， YATES J G. Experimental observations of fluidized beds at elevated pressures［J］. Chemical engineering communications， 1986， 41：133-149.

［8］" 高智雪， 郝振華， 宋素芳， 等. 加壓二維鼓泡床氣固流動特性的數(shù)值模擬［J］. 煤炭轉(zhuǎn)化， 2018， 41（2）：44-50，66.

GAO Zhixue， HAO Zhenhua， SONG Sufang， et al. Numerical simulation of gas-solid flow characteristics in two-dimensional bubble bed under pressure ［J］. Coal conversion， 2018， 41（2）：44-50，66. （in Chinese）

［9］" SNIDER D， BANERJEE S. Heterogeneous gas chemistry in the CPFD Eulerian-Lagrangian numerical scheme （ozone decomposition）［J］. Powder technology， 2010， 199（1）：100-106.

［10］" 羅占. 冷態(tài)雙流化床氣固流動電容層析成像測量及CPFD數(shù)值模擬研究［D］.北京：中國科學(xué)院工程熱物理研究所， 2020.

［11］" 王帥. 流化床內(nèi)稠密氣固兩相反應(yīng)流的歐拉-拉格朗日數(shù)值模擬研究［D］.杭州：浙江大學(xué)， 2019.

［12］" 王洪遠. 基于動態(tài)雙網(wǎng)格的鼓泡流化床傳熱特性的數(shù)值模擬［D］.北京：華北電力大學(xué)， 2021.

［13］" 白玲， 韓晨， 徐云峰， 等. 基于不同曳力模型的鼓泡流化床CFD-DEM數(shù)值模擬與試驗研究［J］. 排灌機械工程學(xué)報， 2022， 40（1）：49-54.

BAI Ling，HAN Chen，XU Yunfeng，et al．Numerical si-mulation and experimental study of CFD－DEM in bubbling fluidized bed basedon different drag models［J］．Journal of drainage and irrigation machinery enginee-ring，2022，40（1）：49-54．（in Chinese）

［14］" MANSOURPOUR Z， KARIMI S， ZARGHAMI R， et al. Insights in hydrodynamics of bubbling fluidized beds at elevated pressure by DEM-CFD approach［J］. Parti-cuology， 2010， 8（5）：407-414.

［15］" LIU D， SONG J， MA J， et al. Gas flow distribution and solid dynamics in a thin rectangular pressurized fluidized bed using CFD-DEM simulation［J］. Powder techno-logy， 2020， 373：369-383.

［16］" 宋加龍. 加壓流化床氣固流動特性研究［D］.南京：東南大學(xué)， 2019.

［17］" GU Y， OZEL A， SUNDARESAN S. A modified cohe-sion model for CFD-DEM simulations of fluidization［J］. Powder technology， 2016， 296：17-28.

［18］" ERGUN S. Fluid flow through packed columns［J］. Chemical engineering progress， 1952， 48：89-94.

（責(zé)任編輯" 黃鑫鑫）