基于CFD-DEM的柱形顆粒固定床反應(yīng)器流場特性數(shù)值模擬

王凱莉,張亞新

新疆大學(xué)化學(xué)化工學(xué)院，新疆 烏魯木齊 830046

固定床反應(yīng)器被廣泛應(yīng)用于各種化工生產(chǎn)過程中，其結(jié)構(gòu)多為圓柱形，內(nèi)部床層由顆粒狀催化劑和催化劑載體填充而成。對于球形顆粒填充床，國內(nèi)外學(xué)者已進(jìn)行深入的研究[1-5]，但在實際生產(chǎn)中，催化劑顆粒形狀多為異形，如圓柱外齒輪形，七筋車輪形和四孔圓柱形等?，F(xiàn)有的研究表明，催化劑顆粒形狀是影響反應(yīng)效果和操作條件的重要因素[6,7]，對固定床床層的流場特性有著極其重要的影響。Dixon[8]研究了顆粒形狀對固定床空隙率的影響，得到球形、空心圓柱體和實心圓柱體的空隙率計算方程。董寶川[9]采用計算流體動力學(xué)（CFD）方法，模擬了固定床內(nèi)部結(jié)構(gòu)對流動性能的影響，認(rèn)為球形比圓柱形顆粒堆積床層更加密實。Nemec等[10]研究了催化劑顆粒形狀對填充床內(nèi)的壓降影響，認(rèn)為Ergun方程只適用于預(yù)測球形顆粒填充床的壓降。Yang等[11]研究了顆粒形狀對結(jié)構(gòu)填充床中的流動和熱傳遞都有重要影響。雖然上述研究揭示了顆粒形狀對床層流場特性的影響，但對于顆粒隨機(jī)堆積固定床的數(shù)值模擬研究，大多局限于球形或顆粒形狀的模型簡化上，而且這些研究沒有涉及柱形催化劑。本工作采用計算流體動力學(xué)-離散單元耦合法（CFD-DEM），通過DEM軟件建立柱狀顆粒實體隨機(jī)堆積床層，而后與Fluent軟件耦合計算，通過床層空隙率、壓力場和速度場分布規(guī)律分析了不同管徑和顆粒直徑的比值（直徑比）對流場特性影響，基于此方法的模擬結(jié)果更加接近生產(chǎn)實際，為反應(yīng)器結(jié)構(gòu)優(yōu)化提供了理論依據(jù)。

1 計算模型

1.1 流體相控制方程

CFD-DEM耦合模型中氣相控制方程主要有A模型、B模型和C模型，具體可參見文獻(xiàn)[12]。本工作采用A模型描述顆粒和流體相的耦合作用，流體相控制方程為：

式中：ε為計算網(wǎng)格內(nèi)的流體空隙率；ρF為流體密度，kg/m3；uF為流體速率，m/s；τ為流體粘性應(yīng)力張量，Pa；g為重力加速度，m/s2；為兩相動量交換原相，即：

1.2 流體-顆粒作用力

采用Free-stream曳力公式計算[13]，即：

式中：ρ為流體密度，kg/m3；A為柱狀顆粒等效體積球的截面積，m2；ν為顆粒與流體之間的相對速度，m/s；cD為阻力系數(shù)，其具體計算公式為：

式中：Re為流體雷諾數(shù)；α為計算流體網(wǎng)格單元的自由體積，m3；L為柱形顆粒等效體積球的接觸直徑，m；η為流體粘度，Pa·s。

2 基于CFD-DEM的床層流場特性數(shù)值模擬

2.1 模型建立

2.1.1 快速填充柱形顆粒模型

以新疆伊犁某企業(yè)實際所用催化劑顆粒為例，統(tǒng)計測量了七孔柱形催化劑顆粒的平均直徑（D）為10.26 mm、平均高度（H）為5.45 mm、平均孔直徑（d）為2.34 mm。運(yùn)用DEM軟件，以球形為基礎(chǔ)單元進(jìn)行填充構(gòu)建柱形顆粒模型。利用workbench平臺先建立如圖1（a）所示的柱形顆粒模型，再通過 fluent網(wǎng)格確定顆粒生成點(diǎn)，借用UDF（User Defined Function）把這些生成點(diǎn)輸出，最后用DEM的顆粒工程API按照已獲得的生成點(diǎn)生成顆粒。此次共填充414顆直徑為1.5 mm的基本球單元顆粒形成柱狀催化劑模型，見圖1（b）。

圖1 柱形催化劑顆粒模擬計算模型Fig.1 Simulation model of cylindrical catalyst particle

2.1.2 固定床反應(yīng)器結(jié)構(gòu)尺寸參數(shù)

圖2為固定床反應(yīng)器結(jié)構(gòu)模型，管長為240 mm，床層高度為150 mm，內(nèi)徑為103 mm。

圖2 固定床反應(yīng)器幾何模型Fig.2 Geometry model of fixed bed reactor

圖3 固定床床層模擬計算模型Fig.3 Simulated computational model of fixed bed reactor

2.1.3 固定床柱狀催化劑床層生成

在完成柱狀催化劑顆粒單元模型的建立后，利用DEM中的Geometry面板中設(shè)置顆粒工廠幾何特征。柱狀催化劑顆粒采用動態(tài)生成的方法，并依據(jù)反應(yīng)器直徑比設(shè)置合適的顆粒總數(shù)和生成速率。柱狀顆粒在自身重力作用下隨機(jī)下落堆積生成顆粒床層，使顆粒填滿整個計算域。建立完成后的顆粒床層模型如圖3所示。

2.2 計算工具條件設(shè)置

2.2.1 DEM設(shè)置

催化劑堆積顆粒采用粗糙毛玻璃材料，固定床的壁面采用碳鋼材料。具體的參數(shù)設(shè)置如下：催化劑顆粒泊松比為 0.25，密度為 2 456 kg/m3，剪切模量為 2.2×108Pa。顆粒碰撞模型采用Hertz-Mindlin(no-slip)模型，顆粒與顆粒恢復(fù)系數(shù)為0.5，顆粒與顆粒靜摩擦系數(shù)為0.154，顆粒與顆粒滾動摩擦系數(shù)為0.1。壁面泊松比為0.3，密度為7 800 kg/m3，剪切模量為7×108Pa，壁面與顆粒碰撞恢復(fù)系數(shù)為0.3，壁面與顆粒靜摩擦系數(shù)為0.154，壁面與顆粒動摩擦系數(shù)為0.1。

2.2.2 Fluent設(shè)置

在Fluent邊界條件設(shè)置中，計算模型選用k-ε湍流模型，模擬為瞬態(tài)，入口條件設(shè)置為速度入口，速率為1 m/s，進(jìn)口與床層之間的距離為30 mm；出口條件設(shè)置為壓力出口，出口壓力為0（表壓），出口與床層的距離為60 mm；采用單求解器、SIMPLEC算法，壁面采用無滑移的標(biāo)準(zhǔn)壁面函數(shù)。模擬過程中選用空氣作為流體相。

2.2.3 耦合模塊設(shè)置

CFD-DEM耦合計算模型有Lagrangian模型和Eulerian模型。本工作屬于氣固兩相流模擬，固相顆粒占用的體積大于總體積的 10%，且顆粒的體積對流體流動的影響不能忽略，因此耦合模型采用Eulerian 模型[14]。

東北四省區(qū)節(jié)水增糧行動項目依靠地下水灌溉的農(nóng)田面積高達(dá)94.5%，水利部對項目建設(shè)合理開發(fā)地下水提出了明確要求。水資源論證是實施好節(jié)水增糧行動項目的重要前期工作，對科學(xué)確定發(fā)展規(guī)模，合理利用有限水資源，緩解水資源供需矛盾，提高水資源利用效率，確保節(jié)水增糧增效有重要意義。本文對項目地下水取水水源、取水影響論證思路及關(guān)鍵技術(shù)問題作了具體介紹。

2.3 模擬結(jié)果數(shù)據(jù)提取及分析

在Fluent后處理中，提取所需的空隙率、壓力和速率數(shù)據(jù)。在固定床軸向（Z方向上）插入直線，提取徑向空隙率和速率分布以及軸向壓力分布，結(jié)果如圖4所示。

圖4 空隙率、速率和壓力模擬結(jié)果Fig.4 Void fraction, velocity and pressure simulation results

由圖 4（a）可知，不同高度的床層空隙率不同，最高值 0.75，最低值為 0.33，但大部分在 0.38和0.53之間。由圖4（b）可知，在同一截面上，速率呈波形變化，在管壁處，其值最小為0，最大值在管壁附近，為3.14 m/s。由圖4（c）可知，壓力值從固定床的進(jìn)口到出口沿著軸向逐步減小，壓降為315.90 Pa。

3 直徑比對流場的影響及分析

在保持床層高度及柱狀催化劑直徑（10.26 mm）不變的情況下，考察不同管徑和顆粒直徑的比值（直徑比，簡寫為N）對流場的影響。

3.1 催化劑床層空隙率分布

圖5為截取了反應(yīng)器軸向剖面（X為0的截面）的流體體積分?jǐn)?shù)云圖。由圖可知，不同直徑比下，床層內(nèi)部均出現(xiàn)局部空隙較大的區(qū)域，這與球形顆粒以自身重力隨機(jī)堆積的床層分布現(xiàn)象相同[15]。當(dāng)N為3時，管壁附近流體空隙體積較大，壁效應(yīng)顯著[16]。當(dāng)N為5時，有空隙局部很小的現(xiàn)象，但顆?？障斗植伎傮w較均勻，空隙率大部分集中在0.39～0.56。當(dāng)N為10時，顆粒分布局部很大或很小的現(xiàn)象較為明顯，流體進(jìn)出口處流體體積分?jǐn)?shù)出現(xiàn)波動，顆粒分布不均勻。當(dāng)N為20時，流體從進(jìn)口流經(jīng)床層80 mm和距出口40 mm時，流體體積分?jǐn)?shù)分層較為明顯，變化趨勢較平穩(wěn)。而在床層80～110 mm，空隙率逐漸變小，顆粒堆積較密實，這使得整個床層顆粒分布均勻性相差較大。因此，N為5時，顆粒分布較為均勻。

圖5 X=0截面軸向顆?？障斗植荚艶ig.5 Distribution contour of axial void fraction in section of X=0

圖6為不同直徑比下，不同床層高度（Z為20，75和130 mm）截面徑向空隙率分布情況。

圖6 不同床層高度的徑向空隙率分布曲線Fig.6 Radial voidage distribution curve at different bed height

由圖可知，同一高度下，不同直徑比下的空隙率分布差別很大。N為3時，空隙率變化波動較大，且空隙率整體值較大；當(dāng)N為5和10時，波動幅度較??；當(dāng)N為20時，徑向空隙率波動較為平緩，但軸向方向上空隙率分層較為明顯。

3.2 床層壓力分布

圖7為不同直徑比下固定床軸向壓力變化情況?？梢钥闯?，壓力從固定床的進(jìn)口到出口沿著軸向方向逐漸減小，而對于不同的直徑比，壓降不同。當(dāng)N為3，5，10和20時，進(jìn)出口的壓降分別為890.96，1 196.30，315.90和91.28 Pa，N為5時，壓降最大，說明N為5時的床層顆粒堆積較密實[17]。

圖7 固定床軸向壓力分布Fig.7 Axial pressure distribution in the fixed bed

3.3 床層速率分布

圖8為不同直徑比下，不同床層高度（Z為20，75和130 mm）截面的速率分布云圖。可以看出，在管壁區(qū)域速率都較小，壁面效應(yīng)較明顯。當(dāng)N為3，10和20時，不同高度下的速率分布云圖變化區(qū)間都較大。當(dāng)N為5時，除了管壁附近，其他區(qū)域速率區(qū)間一致，表明分布較均勻。

圖8 不同直徑比下軸向速率分布云Fig.8 Axial velocity distribution in the fixed bed at different N

圖9為不同直徑比下徑向速率變化規(guī)律。由圖可知，當(dāng)N為3，10和20時，徑向最大速率均在管壁附近，分別為1.47，3.14和18.74 m/s。而N為5時的最大速率在床層中心處附近，最大值為1.63 m/s。分析圖中整體速率變化趨勢可看出，當(dāng)N為3時，速率軸向方向上變化趨勢規(guī)律不明顯。當(dāng)N為5時，速率變化趨勢較平穩(wěn)，中心附近徑向波動多呈直線型。當(dāng)N為10和20時，速率分布由增高到下降再升高，出現(xiàn)波峰，流場分布不均勻，流體流經(jīng)的通道更加曲折。說明N為5的速率流場分布均勻，有利于流體傳遞。

圖9 不同直徑比下徑向速率變化規(guī)律Fig.9 Radial velocity distribution in the fixed bed at different N

4 結(jié) 論

a）基于CFD-DEM耦合法，通過對模擬結(jié)果分析，得出此方法能夠較好地模擬柱形催化劑床層堆積過程及流場分布特性。

b）通過分析不同直徑比下X為0處軸截面的流體體積分?jǐn)?shù)分布云圖和空隙率徑向分布曲線圖，可知不同的直徑比對柱形顆粒隨機(jī)堆積時的影響不同，當(dāng)直徑比為5時，床層顆粒分布較為均勻。

c）通過分析不同直徑比下床層的速度場，可知流體在直徑比為5的床層下流動較平穩(wěn)均勻，波動幅度較小，有利于流體傳遞。當(dāng)直徑比為3時，流體波動不大，但壁效應(yīng)較為顯著，不利于反應(yīng)。

[1]Zhou Xiaoming, Duan Yanjun, Huai Xiulan. 3D CFD modeling of acetone hydrogenation in fixed bed reactor with spherical particles[J].Particuology, 2013, 11(6): 715-722.

[2]Freund H, Zeiser T, Huber F, et al. Numerical simulations of single phase reacting flows in randomly packed fixed-bed reactors and experimental validation[J]. Chemical Engineering Science, 2003, 58(3/6): 903-910.

[3]Guo X. Numerical simulation of flow and heat transfer in a random packed bed[J]. Particuology, 2010, 8(3): 293-299.

[4]Das S, Deen N G, Kuipers J A M. A DNS study of flow and heat transfer through slender fixed-bed reactors randomly packed with spherical particles[J]. Chemical Engineering Science, 2017, 160: 1-19.

[5]Eppinger T, Seidler K, Kraume M. DEM-CFD simulations of fixed bed reactors with small tube to particle diameter ratios[J]. Chemical Engineering Journal, 2011, 166(1): 324-331.

[6]華衛(wèi)琦, 陳豐秋, 呂德偉. 催化劑顆粒的形狀優(yōu)化Ⅰ基本幾何形狀顆粒[J]. 化學(xué)反應(yīng)工程與工藝, 2000, 16(2): 97-102.Hua Weiqi, Chen Fengqiu, Lv Dewei. Shape optimization of catalyst particles I: basic geometry particles[J]. Chemical Reaction Engineering and Technology, 2000, 16(2): 97-102.

[7]華衛(wèi)琦, 陳豐秋, 呂德偉. 催化劑顆粒的形狀優(yōu)化Ⅱ有限長圓柱和常見動力學(xué)[J]. 化學(xué)反應(yīng)工程與工藝, 2000, 16(2): 103-108.Hua Weiqi, Chen Fengqiu, Lv Dewei. Shape optimization of catalyst particles II. Finite length cylinder and common dynamics[J].Chemical Reaction Engineering and Technology, 2000,16(2): 103-108.

[8]Dixon A G. Correlations for wall and particle shape effects on fixed bed bulk voidage[J]. Canadian Journal of Chemical Engineering,1988, 66(5): 705-708.

[9]董寶川. 固定床反應(yīng)器結(jié)構(gòu)及顆粒特性對流動性能的影響[D]. 常州: 常州大學(xué), 2015.

[10]Nemec D, Levec J. Flow through packed bed reactors, part 1: single-phase flow[J]. Chemical Engineering Science, 2005, 60(24):6947-6957.

[11]Yang J, Wang Q, Zeng M, et al. Computational study of forced convective heat transfer in structured packed beds with spherical or ellipsoidal particles[J]. Chemical Engineering Science, 2010, 65(2): 726-738.

[12]任立波. 稠密顆粒兩相流的CFD-DEM耦合并行算法及數(shù)值模擬[D]. 山東大學(xué), 2015.

[13]Hobbs A. Simulation of an aggregate dryer using coupled CFD and DEM methods[J]. International Journal of Computational Fluid Dynamics, 2009, 23(2): 199-207.

[14]王國強(qiáng), 郝萬軍, 王繼新. 離散單元法及其在EDEM上的實踐[M]. 西安: 西北工業(yè)大學(xué)出版社, 2010: 97.

[15]董寶川, 劉雪東, 蘇世卿. 基于CFD-DEM耦合的固定床管式反應(yīng)器流體流動特性數(shù)值模擬[J]. 中國粉體技術(shù), 2015, 21(4): 11-15.Dong Baochuan, Liu Xuedong, Su Shiqing. Numerical simulation of fluid flow characteristics of fixed bed tube reactor based on CFD-DEM coupling[J]. China Powder Science and Technology, 2015, 21(4): 11-15.

[16]Mueller G E. Angular void fraction distributions in randomly packed fixed beds of uniformly sized spheres in cylindrical containers[J].Powder Technology, 1993, 77(3): 313-319.

[17]趙慶國, 廖 暉, 李紹芬. 氣體的溫度和壓力及顆粒形狀對固定床壓降的影響[J]. 化學(xué)反應(yīng)工程與工藝, 2000, 16(1): 1-6.Zhao Qingguo, Liao Hui, Li Shaofen. Effect of temperature and pressure of gas and particle shape on pressure drop of fixed bed[J].Chemical Reaction Engineering and Technology, 2000, 16(1): 1-6.