雙孔射流流化床內顆?；旌咸匦缘碾x散單元法數值模擬

張俊強， 紀 律， 李 斌， 陳海生， 邊禹銘

(1.華北電力大學 能源動力與機械工程學院， 河北保定 071003；2.中國科學院工程熱物理研究所， 北京 100190)

雙孔射流流化床內顆?；旌咸匦缘碾x散單元法數值模擬

張俊強1,2， 紀 律2， 李 斌1， 陳海生2， 邊禹銘1

(1.華北電力大學 能源動力與機械工程學院， 河北保定 071003；2.中國科學院工程熱物理研究所， 北京 100190)

將計算流體力學與離散單元法相結合，采用Fortran語言編程，根據牛頓第三定律實現氣固耦合，從微觀角度剖析了顆粒在流化床內的運動機制.利用Lacey混合指數對流化床內不同特性區域的顆?；旌铣潭冗M行了定量分析，并研究了顆粒混合特性的影響因素，得到顆粒軸向和徑向混合序列圖、氣體和顆粒速度分布以及整床和三區的顆粒混合指數分布.結果表明：顆粒受到空氣射流作用后，隨著混合和偏析的不斷進行，最后達到隨機完全混合狀態；相同表觀氣速下顆粒軸向混合進程快于徑向混合進程；表觀氣速增大，噴泉區尺度增大，使顆粒徑向混合能力得到改善，但對軸向混合能力的影響微弱.

流化床； 軸向混合； 徑向混合； 離散單元法； 混合指數

流化床憑借其特殊的氣固流動特性和優良的結構形式，已在各個工業領域得到了廣泛的應用.流化床內發生強烈的氣固流動及相間反應，因此可作為高溫蓄熱設備，具有換熱效果好和溫度分布均勻等優點.

隨著我國電力行業供給側過剩和可再生能源的迅猛發展，發展電力儲能系統已成為大規模利用可再生能源的有效措施，有助于提高我國常規電力系統的效率、安全性和經濟性[1].流化床以其強大的傳熱和蓄熱能力在儲能領域的探索性應用上受到越來越廣泛的關注.但由于流化床內兩相間作用復雜，人們對流化床的微觀反應過程尚未充分了解.流化床內氣固兩相運動對反應器中傳熱傳質過程具有關鍵作用，而顆?；旌线\動規律對其流動和傳熱過程也有重要影響.

一些學者[2-3]通過實驗手段對流化床內氣固流動特性進行分析，但因操作條件及設備實現難度較大，目前較難獲得微觀層面上的顆粒運動信息.近年來，隨著數值模擬發展不斷加快，諸多學者[4-5]通過數值模擬與實驗研究相結合的方法來深入分析流化床內氣固流動特性.其中離散單元法(DEM)具有獲取微觀信息豐富、模擬精確度高等優勢，在氣固兩相流的模擬研究中得到了越來越多的應用[6-9].

目前，對流化床內顆?；旌蠙C理進行剖析是流化床的研究熱點之一[10-11].Sharma等[12]對鼓泡流化床內生物質與生物碳顆粒的混合與分離行為進行了數值模擬.Olaofe等[13]采用計算流體力學與離散單元法相結合的方法(CFD-DEM)研究了流化床內氣固流動特性,并分析了顆粒的混合和偏析過程.朱潤孺等[14]采用離散單元法研究了鼓泡床內顆粒的軸向和徑向混合特性.但對數值模擬與實驗結果對比，且從微觀層次進行的較為全面的定量研究還較少.

筆者將計算流體力學與離散單元法分別應用于氣相與固相的數值計算，采用Fortran語言自行開發數值模擬程序，對流化床內顆粒的軸向、徑向混合過程進行模擬，并利用實驗結果對其進行驗證，分析了氣體和顆粒在流化床內的流動特性，研究了顆粒的軸向、徑向混合特性及其影響因素.同時采用Lacey混合指數對雙孔射流流化床內顆粒的混合程度進行定量分析，揭示了顆粒由完全分離狀態到隨機完全混合狀態的過程機理.

1 數學模型

1.1 氣固兩相運動模型

流化床內顆粒主要受氣體的曳力、自身重力以及顆粒與顆粒(壁面)之間碰撞力的作用[10].氣相模型采用Navier-Stocks方程，湍流運動采用k-ε兩方程數學模型[10].

1.2 氣固兩相之間的耦合作用

1.2.1 氣相對顆粒的曳力

采用文獻[9]中的曳力模型，氣相對顆粒的曳力如下：

當εg>0.8時，利用Wen&Yu等方程可得

(1)

當εg≤0.8時，采用Ergun等式可得

(2)

(3)

式中：n為空隙率修正因子，一般n的取值為4.65；Cd為單顆粒的曳力系數.

當顆粒雷諾數Rep≤1 000時，

(4)

當Rep>1 000時，

(5)

其中，顆粒雷諾數為

(6)

空隙率為

(7)

式中：ΔV為劃分區域的網格體積；kc為所劃分網格內顆粒的個數.

1.2.2 固相對氣相的反作用力

(8)

式中：Fp為固相對氣相的反作用力.

2 模擬對象及參數

2.1 模擬對象

模擬的床體為長×寬×高=150 mm×4 mm×900 mm的雙噴口矩形截面準三維流化床，床身底部32.5～42.5 mm、107.5～127.5 mm位置分別設置一個空氣進口，空氣進口寬度為10 mm，如圖1所示.

圖1 模擬床體示意圖

床體頂部出口位置采用局部單向化處理.氣相和固相在近壁面處分別設置為無滑移邊界條件和滑移邊界條件.模擬中固相采用2 400個直徑均為4 mm的球形顆粒.

2.2 模擬參數

對不同表觀氣速下雙孔射流流化床內2種組分固相的動態混合過程進行模擬，并分析其影響因素，其中氣相和顆粒的基本參數見表1.

3 模擬結果

3.1 流化床內顆粒的流化過程圖

首先將起始床層的顆粒根據軸向混合、徑向混合均分成2組，以確保顆粒最初為完全分離狀態.圖2和圖3為表觀氣速v=2.0 m/s、計算時間t=0～1.8 s時顆粒的混合序列圖.結合圖2和圖3可知，t=0 s時2種顆粒處于彼此分離狀態，隨著噴動氣流向上流動，顆?；旌喜粩嗌钊?噴動氣流從環隙區滲入噴射區，并將其底部顆粒夾帶到噴泉區，顆粒上行穿過床層表面后分別向兩側擴散.在擴散過程中當顆粒自身重力大于曳力時，顆粒開始下行，并與其他顆粒不斷碰撞混合，緩慢下移至環隙區底部后，又被卷吸至噴射區，重新進入下一循環.顆粒在此循環過程中與周圍顆粒不斷進行混合與擴散，混合程度不斷加深，最后達到顆?；旌系膭討B平衡.由此得出，由氣體攜帶作用引起的顆粒循環對其混合效果有關鍵影響.對比圖2與圖3可知，2種組分顆粒在未達到完全混合前，相同時刻下顆粒軸向混合比徑向混合更充分.以上模擬結果與文獻[10]和文獻[11]中結果一致.

圖2 顆粒軸向混合序列圖

圖3 顆粒徑向混合序列圖

為了驗證以上模型及模擬的準確性，進行了相關實驗研究，測量了相同條件下流化床內顆粒流場瞬時圖.圖4給出了顆粒流動模擬與實驗中瞬時圖的對比.由圖4可以看出，模擬與實驗中床層顆粒軸向、徑向混合時，顆粒由完全分離狀態到混合狀態均具有較好的一致性.

圖4 顆粒流動模擬與實驗中瞬時圖的對比

3.2 氣體速度分布

圖5(a)和圖5(b)分別給出了流化床內不同床高h處的氣體水平速度和垂直速度分布.由圖5(a)可以看出，在床層底部，每個噴口左右兩側的氣體水平速度符號均相反，這是因為進口射流的卷吸作用使氣體從兩側向噴口中心聚集.在左側噴口較高床高位置處，氣體水平速度平均值為正值；在右側噴口較高床高位置處，氣體水平速度平均值為負值，這主要是由于雙噴口位置相對于床體中心更靠近壁面，使得從噴口噴出的氣體在較高床高位置處向壁面方向的運動受限，且在與壁面作用后反向運動.由圖5(b)可以看出，相同床高、靠近2個噴口豎直方向軸線處的氣體速度大于兩側區域的氣體速度，且由于2個噴口處的氣體水平速度存在相互擾動，所以在較高床高位置處，氣體垂直速度最大值向兩側偏移.兩側氣體速度絕對值隨著床高的增加而減小，而中間氣體速度絕對值隨著床高的增加呈增大趨勢.對比圖5(a)與圖5(b)可知，相同床高處氣體垂直速度明顯大于水平速度，氣相對顆粒的曳力作用是顆粒運動的主要影響因素，進而使顆粒的軸向混合比徑向混合更劇烈.

3.3 顆粒速度分布

圖6給出了0～1.8 s內不同床高處的顆粒時均速度分布.由圖6(a)可以看出，由于2個噴口之間存在相互擾動，在床層底部，左側噴口處的顆粒水平速度為負值，此時顆粒往左側移動，右側噴口處的顆粒水平速度為正值，此時顆粒往右側移動；而在較高床高位置處，流動規律呈相反趨勢.由圖6(b)可以看出，2個噴口處的顆粒垂直速度為正值，此時顆粒向床層表面移動.由于2個噴口之間存在較強的橫向擾動，隨著床高的增加，2股主噴射氣流向兩側壁面偏移，在靠近壁面處有連續的氣體作用于顆粒，此時顆粒在主噴射氣流的攜帶作用下向上運動，顆粒上行至床層表面后向中心匯集.中心處的顆粒垂直速度為負值，這是因為中心處氣體速度較小，顆粒自身重力大于曳力，因此顆粒向下運動，運動到床體下部后往2個噴口處匯集并進入各自噴射區，從而構成一個完整的顆粒內循環.對比圖6(a)與圖6(b)可知，相同床高處顆粒垂直速度明顯大于水平速度，從而導致顆粒的軸向混合明顯強于徑向混合.

(a)氣體水平速度

(b)氣體垂直速度

(a)顆粒水平速度

(b)顆粒垂直速度

圖7為顆粒速度矢量圖.由圖7能夠更直觀地得出氣相與固相作用后顆粒在流化床內的運動情況，由該圖得到的顆粒移動進程與由顆粒混合序列圖和顆粒速度分布圖得出的顆粒運動相符，同時也證實了前文中氣固流動特性分析的正確性.

圖7 顆粒速度矢量圖

3.4 混合程度評價

隨著顆?；旌线\動不斷進行，單從直觀上對其混合程度進行分析難以得出可靠的結論.引入Lacey提出的混合指數M[4]，對顆粒的混合程度進行深入研究.選定圖2和圖3中的黑色顆粒為示蹤顆粒，同時將床體平均分成若干個取樣單元，對于有限取樣過程，示蹤顆粒質量濃度標準偏差為

(9)

(10)

混合指數M為

(11)

M為無量綱量，M=0表示2種組分顆粒處于完全分離狀態，M=1表示2種組分顆粒處于隨機完全混合狀態，0

為了更全面地認識流化床內顆?；旌咸匦约捌溆绊懸蛩?，將流化床分成5×30個取樣單元，數據結果輸出步長為0.01 s.圖8給出了表觀氣速為2.0 m/s時顆粒軸向、徑向混合指數的變化.由圖8可以看出，在運動初始階段，混合指數增速較快時顆粒彼此混合迅速；快速混合進行一定時間后，混合指數增速變緩，并在接近1處波動，此時可認為顆?；旌馅呌陔S機完全混合狀態，而波動表明混合過程經歷著混合與偏析，由此達到動態平衡.由圖8還可以看出，顆粒軸向、徑向混合指數變化趨勢區別較大，軸向混合指數起始增速快，在4 s左右達到隨機完全混合狀態；徑向混合指數起始增速慢，在7 s左右達到隨機完全混合狀態.相同時刻下軸向混合進程快于徑向混合進程[11]，這與圖2、圖3中顆粒的混合進程相符.

圖8 顆粒軸向、徑向混合指數M隨時間的變化

通過計算網格來劃分三區，即噴射區、噴泉區和環隙區，每個區域由9個計算網格進行統計計算.圖9給出了該模擬條件下床體軸向、徑向三區混合指數的變化.結合圖8與圖9可知，整床達到隨機完全混合狀態所需時間與環隙區混合指數達到平衡所需時間接近.對比圖9(a)與圖9(b)可知，隨著混合時間的推移，顆粒軸向、徑向混合不同區域的混合指數的變化趨勢相似，但在環隙區存在較大差異.軸向混合時環隙區混合指數達到平衡所需時間較徑向混合時短，這就是整床顆粒軸向混合進程快于徑向混合進程的原因.綜上可知，環隙區內顆粒的混合速度是決定流化床內顆粒達到隨機完全混合狀態的關鍵因素，這與文獻[14]中研究結果相同.對比圖9(a)與圖9(b)還可知，噴泉區混合指數與其他2個區域的混合指數接近，表明噴泉區內顆?；旌陷^為充分，這一結果有別于單噴口流化床[14-15].由于雙孔射流的存在增強了氣相對顆粒的攜帶能力，混合顆粒到達床層頂部后的擴散能力增強，因而噴泉區內顆粒的混合更充分.

(a)軸向混合

(b)徑向混合

圖10給出了表觀氣速分別為2.0 m/s和2.33 m/s時顆?；旌现笖档淖兓?對比圖10(a)與圖10(b)可知，增大表觀氣速對顆粒軸向混合進程的影響較小，但卻有效地加快了顆粒徑向混合進程.表觀氣速增大使顆粒軸向混合速度隨之增大，而顆粒運動形成的床層尺度也相應增大，完成整個顆粒內循環的時間變化較小，導致軸向混合狀態變化較小.而表觀氣速的增大可以有效增大噴泉區尺度，顆粒穿透床層后向兩側的運動速度加快，使顆粒徑向擴散能力增強，即徑向混合能力得到改善.

4 結 論

(1)流化床內相同床高處氣體垂直速度明顯大于水平速度，氣相對顆粒的曳力作用是顆粒流動的主要影響因素，顆粒的軸向混合比徑向混合更劇烈.

(2)噴動氣流將從環隙區滲入噴射區內的顆粒夾帶到噴泉區，顆粒上行到達床層頂部后分別向兩側進行擴散，擴散區域的氣體速度減小，當顆粒所受曳力小于自身重力時，顆粒向上運動速度不斷減小，之后顆粒開始向下運動，顆粒下行期間與周圍顆粒不斷碰撞混合，緩慢下移至環隙區底部后，又被卷吸至噴射區，再次進入下一循環，如此往復，混合程度不斷加強，最后達到顆?；旌蟿討B平衡.

(a)軸向混合

(b)徑向混合

(3)相同表觀氣速下，顆粒軸向混合時達到隨機完全混合狀態所需時間較徑向混合時短.

(4)表觀氣速增大，流化床內噴泉區尺度增大，使顆粒徑向混合能力得到改善，但對軸向混合能力的影響微弱.

[1] CHEN H H, CONG T N, YANG W, et al. Progress in electrical energy storage system: a critical review[J]. Progress in Natural Science, 2009, 19(3): 291-312.

[2] SAIDI M, TABRIZI H B, GRACE J R, et al. Hydrodynamic investigation of gas-solid flow in rectangular spout-fluid bed using CFD-DEM modeling[J]. Powder Technology, 2015, 284: 355-364.

[3] 邵應娟, 金保昇, 鐘文琪, 等. 流化床異重非球顆?；旌咸匦栽囼炑芯縖J]. 工程熱物理學報, 2013, 34(4): 672-675.

SHAO Yingjuan, JIN Baosheng, ZHONG Wenqi, et al. Experimental investigation of irregular particles mixing behavior in fludized bed[J]. Journal of Engineering Thermolphysics, 2013, 34(4): 672-675.

[4] LACEY P M C. Developments in the theory of particle mixing[J]. Journal of Chemical Technology and Biotechnology, 1954, 4(5): 257-268.

[5] 胡顥, 沙春發, 邵應娟, 等. 床料對固廢流化床顆?；旌咸匦缘挠绊慬J]. 東南大學學報(自然科學版), 2013, 43(3): 525-530.

HU Hao, SHA Chunfa, SHAO Yingjuan, et al. Effect of bed material on particle mixing in solid waste fluidized bed[J]. Journal of Southeast University (Natural Science Edition), 2013, 43(3): 525-530.

[6] 任立波, 韓吉田, 趙紅霞. 單沉浸管流化床內離散顆粒數值模擬[J]. 浙江大學學報(工學版), 2015, 49(1): 150-156.

REN Libo, HAN Jitian, ZHAO Hongxia. Numerical simulation of discrete particles in fluidized bed with immersed tube[J]. Journal of Zhejiang University (Engineering Science), 2015, 49(1): 150-156.

[7] 李斌, 紀律. 流化床爐內顆?；旌系碾x散單元法數值模擬[J]. 中國電機工程學報, 2012, 32(20): 42-48.

LI Bin, JI Lü. Numerical simulation of particle mixing in circulating fluidized bed with discrete element method[J]. Proceedings of the CSEE, 2012, 32(20): 42-48.

[8] 紀律. 循環流化床流動及磨損特性的DEM數值模擬[D]. 保定: 華北電力大學, 2011.

[9] LI Tingwen, ZHANG Yongmin, HERNNDEZ-JIMéNEZ F. Investigation of particle-wall interaction in a pseudo-2D fluidized bed using CFD-DEM simulations[J]. Particuology, 2016, 25: 10-22.

[10] 田鳳國, 章明川, 齊永鋒, 等. 流化床軸徑向混合特性的數值研究[J]. 中國電機工程學報, 2006, 26(21): 119-124.

TIAN Fengguo, ZHANG Mingchuan, QI Yongfeng, et al. A numerical investigation on axial/lateral mixing in fluidized beds[J]. Proceedings of the CSEE, 2006, 26(21): 119-124.

[11] 李斌, 宋小龍. 循環流化床內顆?；旌咸匦缘臄抵的M[J]. 動力工程學報, 2013, 33(10): 759-764.

LI Bin, SONG Xiaolong. Numerical simulation on mixing characteristics of particles in circulating fluidized bed[J]. Journal of Chinese Society of Power Engineering, 2013, 33(10): 759-764.

[12] SHARMA A, WANG S B, PAREEK V, et al. CFD modeling of mixing/segregation behavior of biomass and biochar particles in a bubbling fluidized bed[J]. Chemical Engineering Science, 2014, 106: 264-274.

[13] OLAOFE O O, PATIL A V, DEEN N G, et al. Simulation of particle mixing and segregation in bidisperse gas fluidized beds[J]. Chemical Engineering Science, 2014, 108: 258-269.

[14] 朱潤孺, 朱衛兵, 邢力超, 等. 矩形噴動床混合特性的三維數值研究[J]. 中國電機工程學報, 2010, 30(17): 12-16.

ZHU Runru, ZHU Weibing, XING Lichao, et al. A three-dimensional numerical investigation on particle mixing characteristics in rectangular spouted beds[J]. Proceedings of the CSEE, 2010, 30(17): 12-16.

[15] ZHANG Yong, JIN Baosheng, ZHONG Wenqi, et al. DEM simulation of particle mixing in flat-bottom spout-fluid bed[J]. Chemical Engineering Research and Design, 2010, 88(5/6): 757-771.

DEM Simulation on Mixing Characteristics of Particles in Double Jets Fluidized Bed

ZHANGJunqiang1,2,JILü2,LIBin1,CHENHaisheng2,BIANYuming1

(1.School of Energy, Power and Mechanical Engineering, North China Electric Power Uninversity,Baodin 071003, Heibei Province, China; 2. Institude of Engineering Themophsics, Chinese Academy of Sciences, Beijing 100190, China)

Combining the discrete element method (DEM) with computational fluid dynamics, and based on Fortran language programming, a gas-solid coupling model was established according to Newton third law to reveal the motion mechanism of particles in the fluidized bed from the micro-level. The mixing degree of particles in different characteristic regions of the bed was quantitatively analyzed using Lacey mixing index, while factors infuencing the mixing characteristics were studied, during which the axial and radial mixing sequence diagram of particles, the gas and particle velocity distribution, as well as the particle mixing index in the entire bed and at different regions of the bed were obtained. Results show that under the action of air jets, the particles get mixed and segregated, and finally achieve dynamic equilibrium. The axial mixing of particles is faster than radial mixing at the same superficial gas velocity. The scale of fountain area increases with rising superficial gas velocity, while the radial diffusion capacity improves, and the axial diffusion capacity lowers.

fluidized bed; axial mixing; radial mixing; discrete element method; mixing index

2016-03-29

國家國際科技合作專項資助項目(2014DFA60600)；國家自然科學優秀青年基金資助項目(51522605)

張俊強(1991-)，男，河北廊坊人，碩士研究生，主要從事強化傳熱、氣固兩相流數值模擬等方面的研究. 李 斌(通信作者)，男，副教授，博士，電話(Tel．)：0312-7522197；E-mail：binli_871@163.com．

1674-7607(2017)02-0091-07

TK224

A 學科分類號：470.20