擬泡-乳曳力模型在大型高溫費托流化床反應(yīng)器中的CFD模擬研究

王 騰，夏梓洪，陳彩霞

(華東理工大學(xué) 資源與環(huán)境工程學(xué)院，上海 200237)

0 引 言

費托合成技術(shù)是煤間接液化的核心部分，是煤清潔利用的關(guān)鍵技術(shù)之一[1]。根據(jù)目標產(chǎn)物的不同及采用的催化劑、反應(yīng)器形式和反應(yīng)條件不同，費托合成又可分為低溫費托合成和高溫費托合成2種工藝[2]。其中高溫費托合成工藝在獲得大量油品的同時，還可以獲得較多的烯烴、含氧化合物等更高附加值的化學(xué)品，是未來煤制油行業(yè)的重要發(fā)展方向之一[3]。氣固鼓泡流化床由于其氣體和固體顆粒之間傳熱性能好，床內(nèi)溫度均勻，在高溫費托合成反應(yīng)中具有顯著的應(yīng)用優(yōu)勢[4-5]。

目前，對氣固鼓泡流化床中流動、傳遞和混合過程的研究主要依靠試驗完成，諸多關(guān)鍵參數(shù)難以直接測量，流場內(nèi)部特性仍缺乏很好的理解。計算流體力學(xué)(CFD)可以得到氣固流場的詳細信息及其隨時間變化的動態(tài)行為，已成為流化床反應(yīng)器開發(fā)、設(shè)計和放大應(yīng)用的重要輔助工具?；陔p流體模型(TFM)的CFD模擬已經(jīng)成功預(yù)測了氣固流化床內(nèi)的流體力學(xué)特性。然而雙流體模型基于局部平均的假設(shè)，認為單位控制體內(nèi)氣固兩相均勻分布，網(wǎng)格必須足夠小才能正確揭示介尺度(氣泡或顆粒團聚物)的所有細節(jié)[6]。對于氣固鼓泡流化床，網(wǎng)格尺寸至少需要2～4倍顆粒直徑[7-8]。如此細密的網(wǎng)格導(dǎo)致工業(yè)級流化床裝置的計算量遠超現(xiàn)有計算能力。未經(jīng)修正的粗網(wǎng)格計算又會導(dǎo)致顯著的計算誤差[9]。因此，應(yīng)用粗網(wǎng)格模擬大型工業(yè)化裝置時需要結(jié)合非均勻曳力模型對曳力進行修正[10]，以考慮網(wǎng)格尺度內(nèi)未被解析的局部非均勻結(jié)構(gòu)(或稱介尺度結(jié)構(gòu))對氣固流場產(chǎn)生的影響。Wang等[10]采用SGS模型，驗證了使用200倍粒徑粗網(wǎng)格可以得到與試驗結(jié)果一致的固含率分布情況。

氣泡作為氣固鼓泡流化床中最典型的非均勻結(jié)構(gòu)特征，對氣固相間曳力具有決定性影響。要建立可靠的非均勻曳力模型，需考慮氣泡尺寸對曳力的影響。然而，現(xiàn)有的非均勻曳力模型均是在雙流體模型的基礎(chǔ)上對傳統(tǒng)氣固曳力進行修正，以氣含率函數(shù)的形式體現(xiàn)局部非均勻性對曳力的影響，不能完全反應(yīng)氣泡尺寸的影響。實際上，在氣固鼓泡流化床中，氣體分別分布于氣泡相和乳化相中，氣體的不同存在形式使氣固相間作用力的作用機制也不同。乳化相內(nèi)氣體和顆粒間通過直接接觸的方式作用；氣泡與乳化相內(nèi)顆粒間則通過氣泡相和乳化相作用的方式間接作用，且受氣泡尺寸的影響。為了描述氣泡相與顆粒相之間的曳力，本文提出一種基于歐拉多相流框架的新思路——擬泡-乳三相曳力模型(Pseudo bubble-emulsion triple-phase drag model，PBTD)。PBTD模型將氣固鼓泡流化床分為乳化相氣體、乳化相顆粒以及氣泡三相，分別建立守恒方程，體現(xiàn)氣泡的非均勻特性對氣固曳力的影響，乳化相內(nèi)的氣固曳力以及氣泡相與乳化相內(nèi)顆粒的曳力分開考慮。

本文采用開發(fā)的非均勻曳力模型對大型工業(yè)化高溫費托流化床反應(yīng)器進行模擬，考察床層膨脹高度、床層截面氣泡和顆粒體積分數(shù)的分布規(guī)律，研究H2和CO的傳質(zhì)反應(yīng)速率以及氣體組分沿床層高度的分布，并將出口氣體組分與試驗結(jié)果對比，探討該模型在模擬大型工業(yè)化流化床反應(yīng)器中的適用性。

1 模擬方法

1.1 控制方程

本文采用擬泡-乳三相模型將流化床系統(tǒng)分為氣泡相、顆粒相和氣相三相，顆粒相和氣相組成擬乳化相。氣泡相(b)、顆粒相(s)和氣相(g)分別采用一套控制方程進行封閉。

質(zhì)量守恒方程：

(1)

動量方程：

(2)

(3)

(4)

式中，q為g、b、s相；α為體積分數(shù)；ρ為密度；u為速度向量；g為重力加速度；?·τ為黏性力；α?p為靜壓力；?ps為固相壓力；Fg、Fb和Fs為相間作用力項。

應(yīng)力-應(yīng)變張量項τq表示為

(5)

式中，μq為q相剪切黏度；λq為q相體積黏度；I為單位張量。

顆粒相體積黏度項λS[11]為

(6)

顆粒相剪切黏度μs包括動力部分、碰撞部分和摩擦部分，即

μs=μs,col+μs,kin+μs,fr，

(7)

顆粒相碰撞黏度[12]為

(8)

顆粒相動力黏度[13]為

(9)

顆粒相摩擦黏度[14]為

(10)

顆粒相壓力ps方程[11]為

ps=αsρsΘs+2ρs(1+ess)αs2g0,ssΘs，

(11)

其中，ds為顆粒直徑；ess為顆粒碰撞彈性恢復(fù)系數(shù)，取0.9；Θs、g0,ss分別為顆粒溫度、徑向分布函數(shù)；I2D為偏應(yīng)力張量中的第二不變式；φ為內(nèi)摩擦角。

徑向分布函數(shù)[15]為

(12)

式中，αs,max為最大顆粒堆積體積分數(shù)。

顆粒溫度方程：

(13)

(14)

能量傳遞項φs計算公式為

φs=-3βgsΘs。

(15)

式中，βgs為氣固相間動量交換系數(shù)。

1.2 湍流模型

采用混合(Mixture)RNGk-ε模型[16-17]，壁面函數(shù)選擇標準壁面函數(shù)。湍流動能和湍流耗散率控制方程如下：

(16)

(17)

其中：

(18)

(19)

(20)

(21)

Gk,m=μt,m(?um+(?um)T)∶?um，

(22)

式中，ρm為混合密度；μm為分子黏性，um為混合速度；μt,m為湍流黏度；Gk,m為湍動能產(chǎn)生項；αq、ρq、μq、uq分別為q相的體積分數(shù)、密度、黏度和速度；Cμ、C1ε、C2ε為湍流模型參數(shù)；σk、σε為湍流動能k和湍流耗散率ε的湍流普朗特數(shù)。

1.3 曳力模型

氣固流化床中，氣固兩相密度相差很大，其他作用力對氣固流場的影響有限，可忽略不計，只考慮曳力的作用。由于采用了三相假設(shè)，原則上會出現(xiàn)3組相間曳力：氣泡-顆粒曳力、氣泡-乳化相氣體曳力以及乳化相內(nèi)氣體-顆粒曳力。本文在模型推導(dǎo)過程中，將顆粒相和氣相看作“擬乳化相”，考慮氣泡與擬乳化相之間的曳力，并最終將該曳力體現(xiàn)在顆粒和氣泡的動量方程中。在計算氣泡-擬乳化相曳力時采用擬乳化相的密度和黏度，以考慮乳化相中氣體的存在對氣泡和乳化相之間作用力的影響。因此式(2)～(4)中相間作用力項可表示為

Fg=βg-s(ug-us)，

(23)

Fb=βb-s(ub-us)，

(24)

Fs=βb-s(ub-us)+βg-s(ug-us)，

(25)

式中，βg-s、βb-s分別為乳化相內(nèi)氣體-顆粒動量交換系數(shù)、氣泡-顆粒動量交換系數(shù)。

在考慮氣泡與乳化相顆粒相間曳力時，乳化相為連續(xù)相，氣泡為離散相。根據(jù)牛頓第二定律，乳化相對氣泡相的曳力等于氣泡相對乳化相的曳力，因此根據(jù)通用曳力公式，在單位控制體積內(nèi)，單個氣泡對擬乳化相曳力的表達式為

(26)

其中，CDb為氣泡群曳力系數(shù)；ρe為擬乳化相密度；Uslip,i為相間表觀滑移速度，Uslip,i=(ub-ue)(1-αb)。由此可得，單位控制體內(nèi)，所有氣泡與擬乳化相之間的曳力為

(27)

根據(jù)A類雙流體模型曳力系數(shù)公式[13]，氣泡-擬乳化相之間的相間動量交換系數(shù)可表示為

(28)

式中，αe為乳化相體積分數(shù)；ue為乳化相速度。

由于在PBTD模型中，氣泡-擬乳化相曳力最終反映在顆粒相的動量方程中。因此將式(28)中的ub替換為us。氣泡相與顆粒相之間的等效相間動量交換系數(shù)為

(29)

采用Ishii和Zuber[18]提出的公式，即

CDb=CDb0(1-αb)-0.5，

(30)

式中，CDb0為單氣泡曳力系數(shù)，由Darton和Harrison[19]關(guān)系式計算所得。

(31)

(32)

其中，Reb為氣泡雷諾數(shù)；db為氣泡直徑；μe為擬乳化黏度。

db采用Mori和Wen[20]提出的經(jīng)驗公式，即

(33)

db0=0.003 76(u0-umf)2,

(34)

式中，dbm為最大氣泡直徑；db0為初始氣泡直徑；h為床層高度；Dt為床層直徑；u0為入口氣速；umf為最小流態(tài)化氣速。

本文模擬所用反應(yīng)器內(nèi)置大量冷卻盤管，內(nèi)構(gòu)件的存在會對氣固流體動力學(xué)特征產(chǎn)生影響，其復(fù)雜程度也高于空塔流化床反應(yīng)器。Ozawa等[21]發(fā)現(xiàn)在帶有豎直列管的氣固流化床中，A類顆粒和B類顆粒床層內(nèi)部氣泡仍符合隨著床層的升高氣泡尺寸逐漸增大，但受內(nèi)構(gòu)件的影響，平均氣泡尺寸均小于Mori-Wen關(guān)系式預(yù)測值。Ozawa等根據(jù)列管的影響，對Mori-Wen氣泡尺寸公式中的最大氣泡尺寸進行修正，提出以下最大氣泡尺寸修正公式：

dbm=0.59[DtHt(u0-umf)]0.4，

(35)

式中，Ht為列管中軸線距離。

擬乳化相密度ρe和黏度μe受氣體和固體顆粒濃度的共同影響，密度通過擬乳化相內(nèi)部氣固兩相體積加權(quán)平均求得，黏度計算采用Ishii和Zuber[16]的經(jīng)驗公式，具體為

ρe=ρs(1-αeg)+ρgαeg，

(36)

μe=μg(1-αes/αs,max)-2.5αs,max，

(37)

其中，as,max為乳化相中最大顆粒堆積密度，取0.63；αeg和αes分別為乳化相中相對氣含率和固含率，計算公式為

(38)

αes=1-αeg。

(39)

擬乳化相內(nèi)的氣固相間作用力采用傳統(tǒng)Gidaspow曳力模型[12]，需要注意的是，方程中的氣固體積分數(shù)采用擬乳化相內(nèi)的相對體積分數(shù)，相間動量交換系數(shù)βg-s表達為

(40)

(41)

CD=0.44(Res>1 000)，

(42)

(43)

1.4 傳質(zhì)模型

基于擬泡-乳三相曳力模型的假設(shè)，相間傳質(zhì)存在于氣泡相與擬乳化氣相之間，該相間質(zhì)量傳遞主要由氣流穿流項和擴散項兩部分組成。Sit和Grace[22]認為在細小顆粒的流化床中，擴散速率是相間質(zhì)量傳遞速率的控制因素，而在粗顆粒流化床中，氣流穿流是相間傳遞速率的決定性因素。Kunii和Levenspiel[23]將氣泡與乳化相之間的質(zhì)量傳遞過程分為2個步驟：首先氣體通過氣泡傳遞到氣泡暈；其次再從氣泡暈傳遞到乳化相中。并分別用2個交換系數(shù)Kbc和Kce的綜合表示氣泡與乳化相之間的交換總系數(shù)Kbe。本文采用Kunii和Levenspiel提出的質(zhì)量傳遞模型描述氣泡與擬乳化相之間的質(zhì)量傳遞并忽略氣泡暈相，認為氣泡與氣泡暈之間的質(zhì)量交換系數(shù)Kbe等于氣泡與擬乳化相之間的質(zhì)量傳遞系數(shù)，表達式為

(44)

其中，Dg為氣體擴散系數(shù)。等號右側(cè)第1項為氣體穿過氣泡而引起的質(zhì)量交換，即穿流貢獻項(through-flow contribution)，第2項為擴散引起的質(zhì)量交換，即擴散貢獻項(diffusion contribution)。

除了氣泡和乳化相之間的傳質(zhì)，乳化相內(nèi)中氣體混合物之間也存在擴散現(xiàn)象。對于多組分氣體中的擴散，組分i在混合物中的質(zhì)量擴散系數(shù)Di,m可表示為

(45)

其中，Xi為組分i的摩爾分數(shù)；Di,j為雙組分體系的分子擴散系數(shù)，其計算公式為

(46)

1.5 反應(yīng)動力學(xué)模型

F-T合成反應(yīng)是一個多種反應(yīng)并存的復(fù)雜反應(yīng)體系，具有眾多的主反應(yīng)和副反應(yīng)，其發(fā)生的概率隨催化劑種類和操作條件的不同而變化[24]。主反應(yīng)主要生成烷烴類和烯烴類產(chǎn)物，副反應(yīng)包括水煤氣轉(zhuǎn)化(Water gas shift，WGS)反應(yīng)以及生成含氧化合物的反應(yīng)。

主反應(yīng)為：

(47)

副反應(yīng)為：

(48)

(49)

關(guān)于F-T合成反應(yīng)的動力學(xué)，學(xué)者進行了深入廣泛的研究。目前最常用的是考慮合成氣消耗速率的氣體消耗動力學(xué)模型[25-29]。本文采用某公司通過試驗測量數(shù)據(jù)擬合的動力學(xué)模型，具體為

(50)

(51)

其中，RFT、RWGS分別為費托反應(yīng)和水煤氣變換反應(yīng)速率系數(shù)，kmol/(m3·s)；kft、Eft、k1、kwgs、Ewgs為常數(shù)，取值分別為：kft=0.003 224 6，Eft=46 239，k1=0.029，kwgs=1.54，Ewgs=7 392。Keq=0.012exp(4 630/T)；P(H2)、P(CO)、P(H2O)分別為H2、CO和H2O的氣相分壓；R為通用氣體常數(shù)，取8.314 5 J/(K·mol)。

F-T合成反應(yīng)除生成烴類產(chǎn)物外還生成相當數(shù)量的含氧有機物，兩者相對選擇性是一定值，其中含氧有機物選擇性XCHO為

(52)

其中，m為試驗測得的摩爾流率；prod為出口氣體；feed為入口氣體。由于高溫費托合成的含氧有機物組分復(fù)雜，模型采用含量最高的乙醇代表含氧有機物。

1.6 反應(yīng)器幾何模型及計算條件

計算幾何模型基于某公司開發(fā)的4 500 t/a中試規(guī)模高溫費托流化床反應(yīng)器裝置，如圖1所示(L為冷卻盤管長度，D為反應(yīng)器直徑，D′為冷卻盤管直徑)。反應(yīng)器高度為24 000 mm，直徑為900 mm。內(nèi)部冷卻盤管外徑為90 mm，高度為7 000 mm。

圖1 高溫費托流化床反應(yīng)器結(jié)構(gòu)示意

采用六面體結(jié)構(gòu)化網(wǎng)格對幾何模型進行網(wǎng)格劃分，網(wǎng)格數(shù)量約為88萬，屬于平均尺寸約為600倍顆粒直徑的粗網(wǎng)格。為了捕捉壁面對流場的影響，網(wǎng)格在管壁附近進行了局部加密處理，如圖2所示。

圖2 中試規(guī)模高溫費托流化床網(wǎng)格

模擬時氣相和氣泡相采用均勻進氣，進氣速度U0=0.495 m/s。由于催化劑所用顆粒為Geldart A類顆粒，因此氣相進口表觀速度為最小鼓泡速度Umb，為0.007 2 m/s，氣泡相進口表觀速度為(U0-Umb)。入口邊界湍流強度為5%，水力直徑采用幾何直徑。氣相采用無滑移(no-slip)的邊界條件，氣泡相壁面采用自由滑移，顆粒相壁面采用部分滑移(partial slip)，其中鏡面系數(shù)(specularity coefficient)設(shè)置為0.05[30]，顆粒碰撞恢復(fù)系數(shù)設(shè)置為0.9[31]。相間作用力模型和氣泡尺寸模型通過用戶自定義函數(shù)(User defined function，UDF)并耦合到Fluent模擬計算。所用時間步長設(shè)置為0.001 s，每個時間步長迭代20步。模擬時間總長為100 s，對50～100 s模擬瞬時結(jié)果進行數(shù)據(jù)采集并進行時均處理。

模擬采用的物性參數(shù)見表1。試驗測得床層內(nèi)部溫度基本均勻，因此本文假設(shè)床內(nèi)溫度均勻，不存在溫度梯度，模擬不涉及床層與冷卻盤管傳熱，只考慮冷卻盤管對氣固流場的影響。

表1 模擬參數(shù)

2 結(jié)果分析與討論

2.1 總體流化狀態(tài)

圖3為采用非均勻曳力模型時模擬的床層膨脹高度情況，可以看出，床層穩(wěn)定時床層高度為初始床層高度的2倍左右。根據(jù)Cloete[32]提出的床層膨脹高度判據(jù)：床層膨脹高度應(yīng)該在固含率0.05的位置，模擬所得床層最終膨脹高度為8.2 m。由于中試裝置處于高溫高壓的惡劣環(huán)境，無法進行試驗測量，本文采用Cai等[33]提出的經(jīng)驗公式(53)計算床層膨脹高度，該經(jīng)驗公式適用于計算高溫高壓條件下氣固流化床的床層膨脹高度。經(jīng)計算，床層膨脹高度約為8.1 m，模擬床層膨脹高度和經(jīng)驗公式預(yù)測值的相對誤差為1.2%，PBTD曳力模型能夠在粗網(wǎng)格條件下較好地預(yù)測床層的膨脹情況。

(53)

圖3 流化床床層固含率分布

2.2 床層時均相含率分布

床層不同高度處的時均氣泡體積分數(shù)和固體顆粒體積分數(shù)分布云圖如圖4、5所示?？芍鋮s盤管管壁附近，氣泡的體積分數(shù)明顯高于其他部分，固含率則相反，管壁附近固含率低于其他部分，說明氣泡趨于沿著冷卻盤管管壁上升，這與Rudisuli等[34-35]的試驗現(xiàn)象一致。在流化床壁面，氣泡體積分數(shù)降低，顆粒體積分數(shù)升高，這主要是由于顆粒到達流化床床層表面后，氣泡離開床層，顆粒沿壁面下降。表明有內(nèi)構(gòu)件存在的流化床仍然具有固體顆粒濃度邊壁濃、中心稀的非均勻分布特點。隨著高度的增加，床層截面氣泡平均體積分數(shù)逐漸增大，顆粒平均體積分數(shù)逐漸下降，符合流化床顆粒濃度下濃上稀的總體非均勻特性。

圖4 不同床層高度處時均氣泡體積分數(shù)分布云圖

圖5 不同床層高度處時均固含率分布云圖

2.3 傳質(zhì)速率及各組分濃度分布

合成氣的傳質(zhì)速率分布云圖和費托反應(yīng)速率分布云圖如圖6、7所示。可知，CO傳質(zhì)速率≈H2傳質(zhì)速率/2，與反應(yīng)方程中的化學(xué)計量比一致。由圖中還可以看到，伴隨著反應(yīng)器高度的增加，反應(yīng)不斷進行，CO和H2的傳質(zhì)速率和F-T反應(yīng)的反應(yīng)速率都呈逐漸減弱的趨勢。

圖6 H2和CO的傳質(zhì)速率分布

圖7 費托合成反應(yīng)速率分布

床層內(nèi)部氣體成分(CO、H2、[—CH2—]、C2H5OH、H2O、CO2)時均質(zhì)量分數(shù)隨流化床床層高度的變化如圖8所示?？芍S著F-T反應(yīng)的進行，CO和H2被大量消耗，其含量隨床層高度的升高迅速降低。反應(yīng)生成產(chǎn)物[—CH2—]、C2H5OH、CO2、H2O隨著床層高度的升高逐漸增加。床層1.6 m附近，組分曲線明顯出現(xiàn)一次明顯波折，隨后氣體組分反應(yīng)速率變緩慢。這是由于1.6 m后床層開始出現(xiàn)冷卻盤管，使氣固流場發(fā)生明顯改變，氣泡沿管壁上升，氣體傳質(zhì)速率降低，因此反應(yīng)速率隨之降低。

圖8 氣體組分濃度沿床層高度的分布

反應(yīng)器出口氣體組分質(zhì)量分數(shù)與試驗測量數(shù)據(jù)對比如圖9所示，可知主要產(chǎn)物氣體組分質(zhì)量分數(shù)的模擬值與試驗值之間比較接近，主要生成產(chǎn)物(H2O、CO2、C2H5OH和[—CH2—])質(zhì)量分數(shù)相對誤差在1.5%～16.0%，表明PBTD模型基本適用于工業(yè)規(guī)模的反應(yīng)器模擬。

圖9 反應(yīng)器出口氣體組分濃度

3 結(jié) 論

1)采用PBTD模型結(jié)合Ozawa等提出的氣泡尺寸模型，在粗網(wǎng)格條件下對反應(yīng)器進行模擬，得到的床層膨脹高度與文獻公式預(yù)測值的相對誤差為1.2%。PBTD模型能夠呈現(xiàn)氣泡趨于沿著內(nèi)構(gòu)件管壁上升的運動規(guī)律，發(fā)現(xiàn)有內(nèi)構(gòu)件的氣固流化床仍然具有顆粒濃度下濃上稀、邊壁濃中心稀的宏觀非均勻特性。

2)CO傳質(zhì)速率約為H2傳質(zhì)速率的1/2，與反應(yīng)方程中的化學(xué)計量數(shù)一致。模擬結(jié)果顯示，隨著反應(yīng)進行，CO和H2傳質(zhì)和反應(yīng)速率逐漸降低。

3)冷卻盤管的存在影響氣固混合效率，從而影響反應(yīng)速率。模擬出口主要產(chǎn)物氣體組分(H2O、CO2、C2H5OH和[—CH2—])質(zhì)量分數(shù)與試驗測量值相對誤差在1.5%～16.0%，表明本文開發(fā)的基于非均勻假設(shè)的PBTD模型適用于模擬工業(yè)規(guī)模的反應(yīng)器，為大型工業(yè)化流化床反應(yīng)器提供了實用的模擬方案。