氣固流化系統(tǒng)多尺度跨流域EMMS建模

1、氣固流化系統(tǒng)多尺度跨流域EMMS建模氣固流化床內(nèi)的氣體、顆粒與壁面之間的相互作用導(dǎo)致該類(lèi)系統(tǒng)通常呈現(xiàn)復(fù)雜的時(shí)空多尺度結(jié)構(gòu)。隨著流化氣速?gòu)牡偷礁叩淖兓?，氣固系統(tǒng)可能會(huì)依次出現(xiàn)鼓泡、湍動(dòng)、快速流態(tài)化以及稀相輸送等跨流域流動(dòng)特征1。在系統(tǒng)內(nèi)部，局部非均勻性表現(xiàn)為氣泡和團(tuán)聚物的動(dòng)態(tài)生成和湮滅，具有時(shí)間相依性；而整體非均勻性則表現(xiàn)為反應(yīng)器在軸向和徑向呈現(xiàn)有序分布的特征，具有空間相依性1。介尺度結(jié)構(gòu)被認(rèn)為是流化床反應(yīng)器內(nèi)多尺度非均勻性產(chǎn)生的根源，但目前人們對(duì)其認(rèn)識(shí)程度仍不能滿(mǎn)足實(shí)際工藝設(shè)計(jì)的需要。深入理解氣固兩相流介尺度動(dòng)態(tài)結(jié)構(gòu)特性是流化床反應(yīng)器定量放大和優(yōu)化設(shè)計(jì)的基礎(chǔ)和前提。針對(duì)氣固流態(tài)化系統(tǒng)，無(wú)論是顆

2、粒解析的直接數(shù)值模擬（DNS）或是流體解析的離散顆粒模型（DPM），均由于計(jì)算量巨大而難以應(yīng)用于工業(yè)尺度反應(yīng)器的模擬和優(yōu)化。因此，基于粗?；╟oarse-grained）的離散顆粒模擬和基于粗網(wǎng)格（coarse-grid）的連續(xù)介質(zhì)模擬仍然是現(xiàn)階段主要的工業(yè)數(shù)值模擬策略。粗網(wǎng)格處理帶來(lái)的顯著問(wèn)題是平均方法的失效，此時(shí)必須考慮非解析的亞網(wǎng)格結(jié)構(gòu)對(duì)相間作用力和固相應(yīng)力等本構(gòu)關(guān)系的貢獻(xiàn)2-3，且相間曳力關(guān)系的合理構(gòu)建是模擬取得成功的關(guān)鍵2-6。目前主要有三種相間曳力建模方法。（1）經(jīng)驗(yàn)參數(shù)修正法5-8。通過(guò)對(duì)比實(shí)驗(yàn)和模擬數(shù)據(jù)，逆推出曳力系數(shù)的修正關(guān)系。例如，根據(jù)實(shí)驗(yàn)獲得的最小流化速度及空隙率對(duì)傳統(tǒng)

3、均勻曳力模型進(jìn)行校正7；或采用恒定的團(tuán)聚物直徑（dcl）替代均勻曳力模型中的單顆粒直徑（dp）來(lái)計(jì)算實(shí)際曳力系數(shù)8等。這類(lèi)方法的主要缺點(diǎn)是缺乏對(duì)物理機(jī)制的分析，因而模型的可拓展性較差。（2）濾波函數(shù)法3,9-15。通過(guò)濾波函數(shù)對(duì)守恒方程各項(xiàng)進(jìn)行卷積計(jì)算，以獲得空間平均化的控制方程。濾波方程產(chǎn)生的余項(xiàng)來(lái)源于局部非均勻結(jié)構(gòu)，需要采用亞網(wǎng)格模型進(jìn)行封閉。以曳力系數(shù)e(x1,x2,)為例，一般可通過(guò)高分辨率的連續(xù)介質(zhì)模擬、離散顆粒模擬或直接數(shù)值模擬獲得的數(shù)據(jù)進(jìn)行關(guān)聯(lián)9-15。根據(jù)研究體系或選擇參量集xi的不同，眾多的e(x1,x2,)關(guān)聯(lián)函數(shù)被相繼提出。雖然它們最終形式差別很大，但一般均是濾波特征尺寸

4、、固相分率、滑移速度等變量的函數(shù)9-17。近期還進(jìn)一步在模型中考慮了材料性質(zhì)、各向異性、壓力梯度等因素的影響18-20。如何選擇合適的標(biāo)記參量（marker）組合是目前研究的熱點(diǎn)之一20-21。與傳統(tǒng)濾波方法不同，基于湍流的建模方法則是在質(zhì)量和動(dòng)量守恒方程的基礎(chǔ)上嚴(yán)格推導(dǎo)本構(gòu)關(guān)系22-23，并通過(guò)對(duì)濾波余項(xiàng)的權(quán)重分析簡(jiǎn)化封閉關(guān)系的表達(dá)。（3）介尺度結(jié)構(gòu)建模法24-31。通過(guò)在亞網(wǎng)格尺度上對(duì)局部非均勻流動(dòng)結(jié)構(gòu)進(jìn)行分解，以氣泡或團(tuán)聚物為研究對(duì)象建立一系列守恒關(guān)系，并在附加約束條件下通過(guò)解析求解或數(shù)值求解方法獲得曳力系數(shù)的顯式表達(dá)。典型的代表是能量最小多尺度（EMMS）曳力模型25-26。它將系統(tǒng)分

5、解為稀密兩相并分別建立質(zhì)量和動(dòng)量守恒方程，進(jìn)一步引入團(tuán)聚物方程和穩(wěn)定性約束條件，實(shí)現(xiàn)模型的封閉和數(shù)值求解1,27。與經(jīng)驗(yàn)方法及濾波模型不同，EMMS模型通過(guò)量化局部非均勻結(jié)構(gòu)來(lái)更直接地表征操作條件或流域變化時(shí)局部介尺度結(jié)構(gòu)的演變及其對(duì)曳力的影響。EMMS模型的提出最初是為了計(jì)算快速流化床的整體非均勻結(jié)構(gòu)。為了與計(jì)算流體動(dòng)力學(xué)（CFD）相結(jié)合，Yang等26在EMMS模型中引入了顆粒加速度并將其應(yīng)用于計(jì)算網(wǎng)格，從而給出了一種基于結(jié)構(gòu)的EMMS曳力系數(shù)表達(dá)式。此后，該思想被研究者廣泛采用，許多改進(jìn)模型被相繼提出以適應(yīng)更廣泛的操作條件28-37。由于原模型中團(tuán)聚物方程在低顆粒通量、強(qiáng)慣性力和負(fù)顆粒速

6、度情況下并不適用（可能導(dǎo)致dcldp或dcl0）2,30,38，因此很多研究嘗試改進(jìn)團(tuán)聚物方程的表達(dá)式以提高EMMS曳力的普適性13,37-42。有研究指出EMMS求解域過(guò)寬會(huì)導(dǎo)致密相固含率偏離真實(shí)值39，采用經(jīng)驗(yàn)關(guān)聯(lián)式直接關(guān)聯(lián)密相空隙率則可以規(guī)避這一問(wèn)題且可以降低計(jì)算復(fù)雜度13,41-42。另一種思路是尋找新的約束方程替代經(jīng)驗(yàn)性的團(tuán)聚物方程，如Wang等29通過(guò)分析虛擬質(zhì)量力的不同表達(dá)關(guān)聯(lián)了稀相和密相顆粒的加速度，Hu等30則通過(guò)量化稀密相之間的固相質(zhì)量傳遞建立了團(tuán)聚動(dòng)態(tài)平衡方程，基于這些方程改進(jìn)的EMMS曳力被證明具有更好的通用性。近年來(lái)，EMMS思想還被拓展到鼓泡床43-49、湍動(dòng)床50

7、-52、下行床51,53等流域的曳力建模。此外，由于傳統(tǒng)的EMMS曳力模型并不滿(mǎn)足伽利略不變性，因此在改變操作條件時(shí)需要重構(gòu)曳力關(guān)聯(lián)式。所以，Hu等54基于EMMS解集中密相空隙率和穩(wěn)定性條件的關(guān)系提出了一種滿(mǎn)足伽利略不變性的簡(jiǎn)化曳力模型。近年來(lái)，機(jī)器學(xué)習(xí)開(kāi)始被應(yīng)用于粗網(wǎng)格本構(gòu)關(guān)系的關(guān)聯(lián)，這些數(shù)據(jù)驅(qū)動(dòng)的本構(gòu)關(guān)系模型可以在不明顯損失精度的前提下更好地考慮如物性和流動(dòng)等參數(shù)對(duì)曳力的影響55-58。盡管EMMS模型有效考慮了局部非均勻結(jié)構(gòu)對(duì)曳力本構(gòu)關(guān)系的影響，但既往研究大多圍繞靜態(tài)的氣泡和團(tuán)聚物建立模型，網(wǎng)格內(nèi)介尺度結(jié)構(gòu)特征在特定操作參數(shù)和局部空隙率情況下為一定值。這樣的處理方式忽略了氣泡和團(tuán)聚物動(dòng)

8、態(tài)演化的歷史信息，很大可能導(dǎo)致計(jì)算網(wǎng)格內(nèi)實(shí)際的介尺度結(jié)構(gòu)信息和EMMS預(yù)測(cè)值不匹配，因此有必要在模型中進(jìn)一步引入時(shí)間變量并考慮氣泡和團(tuán)聚物的時(shí)空動(dòng)態(tài)演化。群平衡模型（PBM）考慮了離散相的聚并破碎動(dòng)力學(xué)對(duì)數(shù)密度函數(shù)的影響，已經(jīng)被廣泛應(yīng)用于多相體系中離散相尺寸分布的時(shí)空演化特征表征59-61。Hu等62-63通過(guò)量化稀密相質(zhì)量交換建立了團(tuán)聚物的連續(xù)生長(zhǎng)模型以及低通量床中的團(tuán)聚物群平衡方程，并在此基礎(chǔ)上完善了連續(xù)介質(zhì)模型、EMMS曳力和群平衡模型的耦合框架，隨后還拓展了PBM和EMMS耦合模型在鼓泡流化床中的應(yīng)用。由于考慮了介尺度結(jié)構(gòu)時(shí)空動(dòng)態(tài)演化行為對(duì)局部氣泡和團(tuán)聚物尺寸分布的影響，所以PBM模型

9、的引入使得EMMS曳力的修正結(jié)果更加可靠。計(jì)算流體動(dòng)力學(xué)模擬方法在被應(yīng)用于工業(yè)反應(yīng)器時(shí)可能會(huì)難以兼顧效率、精度和成本，但通過(guò)優(yōu)化模擬初場(chǎng)、從特定的整體分布開(kāi)始迭代可以改善上述問(wèn)題64-65。一方面，相對(duì)顆粒均勻分布和緊密堆積的初始狀態(tài)，通過(guò)宏尺度模型對(duì)擬穩(wěn)態(tài)流場(chǎng)進(jìn)行預(yù)測(cè)無(wú)疑更符合實(shí)際工況30。另一方面，人們有時(shí)更關(guān)注反應(yīng)器在不同操作條件下的整體表現(xiàn)，采用大批量的數(shù)值模擬成本高昂且效率低下，此時(shí)穩(wěn)態(tài)建模也可以為反應(yīng)器的設(shè)計(jì)優(yōu)化提供直接依據(jù)66-70。然而，氣固流化床操作方式多變，且伴隨著狀態(tài)多值性和流域轉(zhuǎn)變等復(fù)雜的非線(xiàn)性跨流域特征。隨操作氣速的變化，反應(yīng)器內(nèi)的氣固流動(dòng)還可能經(jīng)歷從移動(dòng)床到快速床等

10、連續(xù)演化過(guò)程。這些特征給流化床反應(yīng)器的宏尺度建模帶來(lái)了挑戰(zhàn)1。經(jīng)驗(yàn)的穩(wěn)態(tài)模型往往很難捕捉流域轉(zhuǎn)變過(guò)程中出現(xiàn)的“噎塞”和“液泛”等跨流域特征，但EMMS思想為氣固系統(tǒng)的跨流域統(tǒng)一建模提供了理論基礎(chǔ)65。Hu等30,67通過(guò)考慮提升管反應(yīng)器內(nèi)顆粒加速效應(yīng)和壁面的影響，在局部和全局穩(wěn)定性條件的約束下將EMMS模型進(jìn)行了軸徑向二維擴(kuò)展。通過(guò)分析不同流域條件下的控制機(jī)制及穩(wěn)定性條件表達(dá)，相繼提出了鼓泡EMMS穩(wěn)態(tài)模型45和并/逆流EMMS穩(wěn)態(tài)模型68-69。在此基礎(chǔ)上，Liu等69-70建立了基于EMMS的全流域統(tǒng)一建模方法并繪制了操作相圖。Tu等71進(jìn)一步指出，當(dāng)對(duì)多回路循環(huán)流化床系統(tǒng)進(jìn)行全回路模擬時(shí)

11、，必須恰當(dāng)?shù)乜紤]各子反應(yīng)器模塊的邊界條件及單元操作之間的強(qiáng)耦合作用。EMMS操作相圖為循環(huán)流化床反應(yīng)器內(nèi)的流域識(shí)別和量化表征提供了可行途徑。通過(guò)考慮各模塊之間的壓降平衡和全床物料守恒關(guān)系，所謂的EMMS全回路建模方法可以實(shí)現(xiàn)對(duì)該類(lèi)反應(yīng)器整體穩(wěn)態(tài)動(dòng)力學(xué)的準(zhǔn)實(shí)時(shí)預(yù)測(cè)70,72。本文將重點(diǎn)介紹和總結(jié)本課題組在氣固系統(tǒng)介尺度建模和EMMS模型擴(kuò)展方面的工作。首先介紹基于團(tuán)聚物或氣泡的EMMS曳力模型的改進(jìn)和發(fā)展；第二部分論述介尺度結(jié)構(gòu)時(shí)空動(dòng)態(tài)演化行為對(duì)本構(gòu)關(guān)系的影響，以及群平衡模型在氣固流態(tài)化系統(tǒng)中的應(yīng)用；第三部分介紹基于EMMS的反應(yīng)器宏尺度建模方法，包括模型的擴(kuò)展、相圖繪制和全回路建模等；最后是對(duì)

12、未來(lái)工作的展望。1 EMMS曳力模型如圖1所示，EMMS模型將非均勻氣固系統(tǒng)分解為顆粒富集的密相和氣體富集的稀相，兩者之間通過(guò)相界面相互作用1。假設(shè)兩相內(nèi)的氣固相互作用可以用均勻模型描述，則特定操作條件系統(tǒng)狀態(tài)可以用8個(gè)參數(shù)X=Upf,Uf,f,Upc,Uc,c,f,dcl來(lái)定義。通過(guò)對(duì)稀密兩相分別建立質(zhì)量和動(dòng)量守恒方程，再加上兩者之間的壓降平衡方程及團(tuán)聚物直徑關(guān)聯(lián)式，總共6個(gè)守恒方程得以建立。Li等1,27認(rèn)為在氣固流態(tài)化系統(tǒng)存在兩種控制機(jī)制：當(dāng)流體控制時(shí)，Wst=min；顆粒控制時(shí)，=min。兩者的競(jìng)爭(zhēng)協(xié)調(diào)關(guān)系構(gòu)成了系統(tǒng)的穩(wěn)定性條件，即minWstg或Nst=fX=min(1)圖1圖1EM

13、MS模型的物理圖譜和表達(dá)式65Fig.1Physical concept and formulation of the EMMS model65以上穩(wěn)定性條件和守恒方程共同構(gòu)成EMMS模型。在給定的宏尺度操作條件下時(shí)，X可以通過(guò)數(shù)值方法求解。然而，由于EMMS模型最初用于描述流化床系統(tǒng)的整體流動(dòng)結(jié)構(gòu)，將其拓展到計(jì)算網(wǎng)格求解有效曳力時(shí)至少需要解決和考慮以下問(wèn)題：（1）局部結(jié)構(gòu)遠(yuǎn)離平衡態(tài)，模型方程如何考慮顆粒加速度的影響；（2）原穩(wěn)定性條件在拓展到其他流域時(shí)應(yīng)該如何改進(jìn)；（3）曳力關(guān)聯(lián)式的形式及其計(jì)算方法。1.1 改進(jìn)的基于團(tuán)聚物的EMMS曳力EMMS模型基于能量平衡假設(shè)建立了團(tuán)聚物方程，但該方程

14、難以適用于低通量、高加速度和順重力情形，這限制了模型在曳力計(jì)算中的應(yīng)用。為了改進(jìn)對(duì)顆粒團(tuán)聚動(dòng)力學(xué)的描述以提高EMMS曳力模型的普適性，Hu等30通過(guò)量化團(tuán)聚物和稀相顆粒之間的質(zhì)量交換關(guān)系，提出了團(tuán)聚物動(dòng)態(tài)平衡方程，并以此替換原團(tuán)聚物方程對(duì)EMMS模型進(jìn)行了改進(jìn)。如圖2(a)所示，通常情況下，流體計(jì)算網(wǎng)格中呈現(xiàn)稀密相共存的非均勻兩相結(jié)構(gòu)。在團(tuán)聚物存續(xù)時(shí)間內(nèi)，稀相和密相顆粒之間的質(zhì)量交換可以認(rèn)為由兩個(gè)可逆的過(guò)程組成。如圖2（b）所示，對(duì)于球形顆粒物，在迎風(fēng)面，稀相顆粒和團(tuán)聚物發(fā)生非彈性碰撞導(dǎo)致顆粒被聚團(tuán)捕捉；而在背風(fēng)面，團(tuán)聚物表面顆粒在隨機(jī)脈動(dòng)作用下進(jìn)入稀相。兩者滿(mǎn)足動(dòng)態(tài)平衡時(shí)可得圖2圖2EMMS

15、曳力模型中亞網(wǎng)格非均勻性（a）和稀密相質(zhì)量傳遞過(guò)程（b）Fig.2Sub-grid heterogeneity (a) and solids mass transfer between the dilute and dense phases (b) in the EMMS drag modelUpf1-f-Upc1-c1-f=41-c3c(2)其中，w是常數(shù)，本研究取1/12。密相顆粒溫度c的計(jì)算可參考文獻(xiàn)30。根據(jù)圖2(a) 可以分別對(duì)稀相和密相建立守恒方程。連續(xù)性方程同原EMMS模型（圖1）一致，但由于局部結(jié)構(gòu)往往遠(yuǎn)離平衡態(tài)，顆粒的加速效應(yīng)常不能忽略，稀密兩相的力平衡方程應(yīng)改寫(xiě)為Ff=34

16、CDf1-fdpgUsf2=1-fp-gg+af(3)fFc+Fi=f34CDc1-cdpgUsc2+34CDifdclgUsi2=f1-cp-gg+ac(4)式中，F(xiàn)f、Fc分別表示單位體積稀相和密相中顆粒所受曳力；Fi表示單位床層中稀相繞流聚團(tuán)產(chǎn)生的界面曳力。三者通過(guò)壓降平衡方程（Ff+Fi/1-f=Fc）相關(guān)聯(lián)。相應(yīng)地，在穩(wěn)定性條件的建立過(guò)程中需考慮顆粒加速度對(duì)懸浮輸送能耗的影響，即Nst=1p1-gFcUcf+FfUf1-f+FiUf1-fmin(5)非均勻性因子定義為EMMS曳力系數(shù)和均勻曳力系數(shù)的比值，即Hd=Ehomo=gFcf+Ff1-f+FihomoUg/g-Up/1-g(6

17、)理論上，在CFD模擬中應(yīng)使用當(dāng)?shù)鼐W(wǎng)格內(nèi)的解析參數(shù)（up，up，g）來(lái)在線(xiàn)求解EMMS模型，這顯然需要消耗巨量的計(jì)算資源。因此，在實(shí)際使用中通常固定氣固速度，在區(qū)間mf, 1.0遍歷空隙率g，通過(guò)數(shù)值方式來(lái)求解EMMS模型，并將Hd關(guān)聯(lián)成為g的函數(shù)或系數(shù)矩陣以與CFD耦合。此處暗含的假設(shè)是滑移速度對(duì)Hd的影響遠(yuǎn)不如空隙率明顯。值得說(shuō)明的是，原EMMS模型中的團(tuán)聚物方程在Up 0時(shí)的求解結(jié)果為負(fù)值，因此無(wú)法適應(yīng)如湍動(dòng)流化床和下行床等的計(jì)算。Hu等51采用上述改進(jìn)EMMS曳力對(duì)固相通量為零的湍動(dòng)流化床進(jìn)行了模擬，仍能獲得滿(mǎn)意的預(yù)測(cè)結(jié)果。對(duì)于氣固并流下行床，盡管氣固皆沿重力方向運(yùn)動(dòng)，且氣固曳力的作用

18、方向在入口處和充分發(fā)展段并不一致，但可以通過(guò)定義符號(hào)因子=ug-up/ug-up建立適用于各發(fā)展段的統(tǒng)一形式的守恒方程組，模型采用穩(wěn)定性條件Nst/NT=min進(jìn)行封閉68。圖3展示了采用改進(jìn)EMMS曳力和Wen-Yu曳力計(jì)算得到的瞬態(tài)顆粒濃度云圖對(duì)比?？梢?jiàn)在不同情況下，改進(jìn)模型均更好地預(yù)測(cè)了反應(yīng)器內(nèi)的氣固非均勻流動(dòng)現(xiàn)象，而均勻曳力由于高估了氣固相間曳力而預(yù)測(cè)結(jié)果不佳。圖3圖3改進(jìn)EMMS曳力（a）,(c),(e),(g)和均勻曳力（b）,(d),(f),(h)對(duì)快速床提升管(a)、湍動(dòng)床（c）、變徑提升管(e)）及下行床(g)的預(yù)測(cè)結(jié)果30,51,73Fig.3Snapshot of sol

19、ids concentration in the fast fluidized riser (a), turbulent bed (c), tapered riser (e) and co-current downer (g) by using the improved EMMS drag model (a),(c),(e),(g) and the homogeneous drag model (b),(d),（f）,(h)30,51,731.2 具有伽利略不變性的EMMS曳力傳統(tǒng)的EMMS曳力計(jì)算依賴(lài)系統(tǒng)宏觀操作參數(shù)（如Ug、Gs等），構(gòu)建滿(mǎn)足伽利略不變性的EMMS曳力模型可以避免當(dāng)操作條件

20、改變時(shí)必須重構(gòu)Hd系數(shù)關(guān)聯(lián)式的問(wèn)題，并可以適應(yīng)具有復(fù)雜幾何結(jié)構(gòu)的流化床反應(yīng)器的曳力計(jì)算。建立滿(mǎn)足伽利略不變性EMMS曳力模型的重要步驟是將方程中的絕對(duì)速度項(xiàng)轉(zhuǎn)化為相對(duì)速度項(xiàng)。根據(jù)滑移速度的定義，式（2）可以重組為Usif1-f-Usff=43c1-c1-f(7)借助相間壓降平衡方程消去式（3）、式（4）中的Fi可得稀密相力平衡方程（k=c代表密相；k=f代表稀相）34CDk1-kdpUsk2=p-ggf1-kac+g+1-f1-faf+g(8)為了獲得滿(mǎn)足伽利略不變性的穩(wěn)定性條件的等價(jià)表達(dá)，Hu等54分析了不同操作條件下EMMS模型的解集，發(fā)現(xiàn)當(dāng)氣固物性不變時(shí)，密相空隙率可以近似表達(dá)為平均空隙

21、率的函數(shù)。對(duì)于典型A類(lèi)顆粒（dp= 60 m,p= 1490 kg/m3）有Nstminc=minmax0.148+0.853g1.56,mf,g(9)質(zhì)量守恒方程（氣相、固相和固含率）、模型基本假設(shè)（fmax,af0）及式（7）式（9）共同構(gòu)成滿(mǎn)足伽利略不變性的EMMS曳力模型?？梢?jiàn)，相對(duì)原來(lái)的非線(xiàn)性?xún)?yōu)化求解，上述模型的求解復(fù)雜度大大下降，可以容易地構(gòu)建Hd關(guān)于（，Re）的關(guān)聯(lián)式或在模擬過(guò)程中直接對(duì)模型進(jìn)行在線(xiàn)求解。圖4(a) 展示了采用該模型計(jì)算得到的非均勻因子，其隨空隙率和Reynolds數(shù)的變化符合物理規(guī)律54。圖4(b)進(jìn)一步比較了本模型、改進(jìn)EMMS模型和實(shí)驗(yàn)數(shù)據(jù)的對(duì)比?？梢?jiàn)式（9

22、）的簡(jiǎn)化并未降低模型的準(zhǔn)確性，但模型的普適性有了很大提高。圖4圖4滿(mǎn)足伽利略不變性的EMMS曳力模型54Fig.4EMMS drag model following Galilean invariance541.3 改進(jìn)的基于氣泡的EMMS曳力EMMS曳力模型也被拓展到氣固鼓泡流態(tài)化系統(tǒng)43-44,46,49-50。但由于經(jīng)驗(yàn)性的氣泡關(guān)聯(lián)式往往基于特定的實(shí)驗(yàn)條件獲得而難以被推廣到其他情況，因而Liu等45提出了采用改進(jìn)穩(wěn)定性條件來(lái)約束氣泡尺寸的方法。如圖5所示，氣固鼓泡流態(tài)化系統(tǒng)可以分解為氣泡相和乳化相，共采用9個(gè)結(jié)構(gòu)參數(shù)X=Ub,fb,b,db,ab,Uge,Upe,e,ape來(lái)定義。但為了

23、簡(jiǎn)化計(jì)算，一般假設(shè)氣泡中不含顆粒，即b=1.0。分別以顆粒和氣泡為研究對(duì)象可以建立一系列連續(xù)性方程和力平衡方程。與上述基于團(tuán)聚物的EMMS模型不同，EMMS鼓泡模型一開(kāi)始就考慮了顆粒和氣泡的加速效應(yīng)，且在穩(wěn)定性條件的表達(dá)上也有很大的變化。這是因?yàn)闅馀菰陂L(zhǎng)大過(guò)程中會(huì)因?yàn)橄嘟缑娴恼龎毫Χa(chǎn)生能量耗散，因此在計(jì)算單位質(zhì)量床層的懸浮輸送能耗（Nst）時(shí)，除了曳力做功外，還應(yīng)額外考慮氣泡克服乳化相壓力膨脹所做的體積功。系統(tǒng)的穩(wěn)定性條件定義為相對(duì)懸浮輸送能耗趨于最小，其具體形式見(jiàn)圖5。改進(jìn)后的穩(wěn)定性條件中包含了氣泡直徑的影響，這為氣泡直徑的取值提供了自然約束，避免了經(jīng)驗(yàn)關(guān)聯(lián)式的引入。圖5圖5EMMS鼓泡曳

24、力模型的物理圖解和公式Fig.5Physical diagram and formulation of the bubbling EMMS drag model在鼓泡床中，非均勻性因子可以表達(dá)為Hd=Ehomo=gFde1-fb+FdihomoUg/g-Up/1-g(10)鼓泡EMMS曳力的計(jì)算及與CFD耦合的方式同1.1節(jié)類(lèi)似，此處不再贅述。如圖6所示，相對(duì)均勻曳力模型，EMMS鼓泡曳力模型可以更好地預(yù)測(cè)反應(yīng)器的非均勻結(jié)構(gòu)，與實(shí)驗(yàn)數(shù)據(jù)呈現(xiàn)更好的一致性，這也說(shuō)明了模型的合理性。關(guān)于模型方程的推導(dǎo)和實(shí)驗(yàn)驗(yàn)證的詳細(xì)信息可以見(jiàn)文獻(xiàn)45。圖6圖6采用EMMS鼓泡曳力和均勻曳力對(duì)A類(lèi)顆粒鼓泡床的模擬結(jié)果

25、對(duì)比45Fig.6Comparison between the bubbling EMMS model and homogeneous model in predicting the gas-solid bubbling beds of Geldart A particles452 介尺度結(jié)構(gòu)時(shí)空動(dòng)態(tài)演化的群平衡建模2.1 基于團(tuán)聚物動(dòng)態(tài)演化的群平衡建模在傳統(tǒng)的EMMS曳力模型中，所取團(tuán)聚物直徑可以認(rèn)為是當(dāng)前操作條件（Ug,Up,）下的最概然值或統(tǒng)計(jì)均值，而團(tuán)聚物的動(dòng)態(tài)特征及尺寸分布（CSD）信息未被考慮。為了定量描述團(tuán)聚物的時(shí)空動(dòng)態(tài)分布及其對(duì)網(wǎng)格內(nèi)曳力的影響，Hu等62提出了基于團(tuán)聚物的群平

26、衡方程，并進(jìn)一步建立了CFD、PBM及EMMS模型的耦合計(jì)算框架。以團(tuán)聚物直徑（L）為內(nèi)坐標(biāo)，團(tuán)聚物尺寸演化的群平衡方程可以表達(dá)為tn(L,x,t)+xuxL,x,tn(L,x,t)+LuLL,x,tn(L,x,t)=BL,x,t+DL,x,t(11)式中，n(L,x,t) 表示尺寸為L(zhǎng)的團(tuán)聚物在t時(shí)刻和x位置的數(shù)密度；B(L,x,t) 和D(L,x,t) 分別表示由于聚并和破碎效應(yīng)導(dǎo)致的團(tuán)聚物生成項(xiàng)和消亡項(xiàng)，一般表達(dá)為BL,x,t=L220LaL3-31/3,L3-32/3nL3-31/3,tn,td+LpgLn,td(12)DL,x,t=nL,t0aL,n,td+gLnL,t(13)由于流

27、化床中顆粒湍流及團(tuán)聚物聚并、破碎的機(jī)制尚不明確，因此初步研究主要聚焦于團(tuán)聚物和稀相顆粒的碰撞和離析等連續(xù)演化過(guò)程。通過(guò)分析單位時(shí)間顆粒和團(tuán)聚物的非彈性碰撞及稀相流體對(duì)聚團(tuán)表面顆粒的侵蝕效應(yīng)，Hu等62建立了團(tuán)聚物的連續(xù)生長(zhǎng)方程，代入PBM可得tnc(L,x,t)+xuxx,tnc(L,x,t)+L(upf-upc)(1-f)2(1-c)-23cnc(L,x,t)=0(14)式(14) 理論上更適用于床層較稀的工況。通過(guò)離散法或矩方法可以對(duì)上述方程進(jìn)行求解，進(jìn)而可以獲得團(tuán)聚物的數(shù)密度函數(shù)及網(wǎng)格內(nèi)的平均直徑分布。前述信息可以用來(lái)封閉EMMS曳力模型中的團(tuán)聚物直徑項(xiàng)。圖7展示了CFD-PBM-EMM

28、S的耦合計(jì)算框架：首先利用CFD模型計(jì)算得到的t時(shí)刻所有網(wǎng)格內(nèi)的相速度、空隙率、湍流特征量并傳遞給PBM；求解PBM可以得到t時(shí)刻網(wǎng)格內(nèi)團(tuán)聚物的尺寸分布（CSD）；將以上CFD和PBM的求解信息傳遞給EMMS曳力模型，進(jìn)而獲得粗網(wǎng)格的有效曳力系數(shù)，并提高CFD的計(jì)算精度。圖7圖7CFD、PBM和EMMS模型的耦合Fig.7Integration of the CFD, PBM and EMMS models耦合模型可以有效避免實(shí)際團(tuán)聚物直徑和傳統(tǒng)EMMS預(yù)測(cè)值不一致的現(xiàn)象，具有更好的自洽性。更重要的是，該算法可以直接計(jì)算反應(yīng)器內(nèi)不同位置處團(tuán)聚物的尺寸分布隨時(shí)間的演化規(guī)律，這對(duì)反應(yīng)器的設(shè)計(jì)和優(yōu)化

29、無(wú)疑具有重要的意義。圖8(a) 展示了基于PBM修正的EMMS曳力校正系數(shù)，可見(jiàn)該模型合理地預(yù)測(cè)了團(tuán)聚物尺寸這一關(guān)鍵參數(shù)對(duì)Hd的影響。通過(guò)圖8(b)可見(jiàn)，相對(duì)均勻曳力和傳統(tǒng)EMMS曳力，本模型可以更好地預(yù)測(cè)流化床內(nèi)的非均勻分布。另外，群平衡模型還可以對(duì)團(tuán)聚物的尺寸分布進(jìn)行直接預(yù)測(cè)，如圖8(c)所示。圖8圖8CFD-PBM-EMMS耦合模型在快速流態(tài)化中的應(yīng)用62Fig.8Implement of the CFD-PBM-EMMS model in the fast fluidized beds622.2 基于氣泡動(dòng)態(tài)演化的群平衡建模在氣固鼓泡流態(tài)化中，氣泡自分布板上部生成之后，往往經(jīng)歷復(fù)雜的聚

30、并、破碎等動(dòng)態(tài)演化過(guò)程74-75。很多研究采用實(shí)驗(yàn)和經(jīng)驗(yàn)建模的方式對(duì)氣泡的尺寸分布進(jìn)行研究和建模，不僅難以適應(yīng)各種復(fù)雜操作條件，也無(wú)法對(duì)氣泡的動(dòng)態(tài)時(shí)空分布進(jìn)行準(zhǔn)確計(jì)算。采用這些經(jīng)驗(yàn)方程進(jìn)行封閉也在很大程度上影響了曳力模型的準(zhǔn)確性和合理性。根據(jù)Hu等63的研究，氣泡的聚并主要由大小氣泡的上升速度差異和氣泡的尾流加速引起的，因此根據(jù)兩者的物理過(guò)程可分別建立其聚并核函數(shù)。大小為di和dj的兩個(gè)氣泡的聚并頻率可以表達(dá)為adi,dj=0.54di+dj2ui-ujr+sin2w/2urel,ww(15)式中，r和w分別表示由于氣泡追趕和尾流效應(yīng)導(dǎo)致的聚并效率63；w表示尾流的平均張角；urel,w表示由

31、于尾流誘導(dǎo)產(chǎn)生的碰撞速度，實(shí)驗(yàn)發(fā)現(xiàn)該附加速度和碰撞發(fā)生時(shí)的氣泡上升速度大小相當(dāng)。氣泡的破碎機(jī)制則更加復(fù)雜，一般認(rèn)為是頂部顆粒快速下落導(dǎo)致的氣泡切割和分裂。通過(guò)對(duì)實(shí)驗(yàn)數(shù)據(jù)回歸發(fā)現(xiàn)氣泡的破碎頻率和其直徑的2.9次冪成正比。因此，直徑為di的氣泡破碎頻率可以通過(guò)式(16)關(guān)聯(lián)獲得gdi=g(db*)gdidb*2.9(16)式中，db*表示氣泡平衡直徑。而破裂導(dǎo)致的子氣泡概率分布函數(shù)則采用實(shí)驗(yàn)統(tǒng)計(jì)的擬合函數(shù)。體積為V的氣泡破碎生成的子氣泡近似滿(mǎn)足以下高斯分布VV=2.13Vexp-V/V-0.50.2712(17)將式（15）式（17）代入式（11）式（13）可以得到氣泡的群平衡方程。采用和圖7一致

32、的建模思路建立鼓泡床的CFD-PBM-EMMS模型。如圖9(a) 所示，借助PBM封閉的鼓泡EMMS模型合理地考慮了氣泡直徑對(duì)非均勻曳力的影響。在氣泡直徑趨于0時(shí)，曳力自然趨于均勻曳力。圖9(b)是模型預(yù)測(cè)結(jié)果和實(shí)驗(yàn)的對(duì)比，可見(jiàn)該模型可以較好地預(yù)測(cè)氣泡尺寸沿床層高度的分布規(guī)律。圖9圖9CFD-PBM-EMMS耦合模型在鼓泡流態(tài)化中的應(yīng)用63Fig.9Implementation of the CFD-PBM-EMMS model in the bubbling fluidized beds633 基于EMMS的全流域操作相圖和全回路穩(wěn)態(tài)建模3.1 穩(wěn)態(tài)EMMS模型的擴(kuò)展和相圖繪制研究表明氣固流

33、態(tài)化和物質(zhì)相變具有相似的多狀態(tài)特征，隨操作條件改變呈現(xiàn)不同的非均勻分布和連續(xù)的流域轉(zhuǎn)變76-78。量化流化床的宏尺度穩(wěn)態(tài)動(dòng)力學(xué)不僅可以快速了解反應(yīng)器的操作狀況，亦可以作為初始分布加速工業(yè)級(jí)反應(yīng)的數(shù)值模擬過(guò)程。Hu等67發(fā)現(xiàn)在提升管底部團(tuán)聚物具有強(qiáng)烈的碰撞、聚并作用，而在頂部團(tuán)聚物組織松散，稀密相差異很小。在此基礎(chǔ)上，通過(guò)引入顆粒加速度建立了EMMS軸向模型，首次實(shí)現(xiàn)了對(duì)提升管反應(yīng)器中軸向S形非均勻分布的EMMS預(yù)測(cè)67。根據(jù)軸向模型獲得各截面處的平均空隙率和氣固速度信息，Hu等30進(jìn)一步改進(jìn)和完善了徑向EMMS模型，并通過(guò)函數(shù)逼近的方法實(shí)現(xiàn)了對(duì)提升管中徑向環(huán)核結(jié)構(gòu)的數(shù)學(xué)預(yù)測(cè)。圖10展示了徑向E

34、MMS模型的物理框架。自上而下，在截面尺度，各物理量必須滿(mǎn)足截面穩(wěn)定性條件和邊界約束的守恒方程；在單元介尺度，各徑向微元體內(nèi)參數(shù)服從于穩(wěn)定性條件和質(zhì)量動(dòng)量守恒方程；而在顆粒微尺度，氣固相互作用在稀密兩相內(nèi)呈現(xiàn)均勻結(jié)構(gòu)，可以用均勻曳力計(jì)算。與提升管反應(yīng)器類(lèi)似，氣固下行床中也存在顆粒團(tuán)聚現(xiàn)象。因此，上述軸徑向建模思路也被應(yīng)用到氣固并/逆流下行床的計(jì)算。將氣固并/逆流下行床進(jìn)行結(jié)構(gòu)分解，分別建立各相的質(zhì)量和動(dòng)量守恒方程以及不同發(fā)展階段的穩(wěn)定性條件，最終可建立EMMS并/逆流下行床數(shù)學(xué)模型68-69。至于鼓泡流態(tài)化系統(tǒng)，氣泡的形成和演化代替了提升管中顆粒團(tuán)聚物的聚并和破碎而成為系統(tǒng)的主要非均勻結(jié)構(gòu)特征

35、。假設(shè)乳化相的顆粒加速度為0，則圖5所示方程可以直接用于鼓泡床穩(wěn)態(tài)動(dòng)力學(xué)的求解。模型無(wú)須經(jīng)驗(yàn)氣泡直徑的輔助就可對(duì)鼓泡床動(dòng)力學(xué)進(jìn)行預(yù)測(cè)45。更多研究細(xì)節(jié)和實(shí)驗(yàn)驗(yàn)證可參考文獻(xiàn)30,45,67-70。圖10圖10提升管徑向非均勻分布的EMMS穩(wěn)態(tài)建模64-65Fig.10Modelling of radial hydrodynamics of risers by using the EMMS model64-65EMMS模型在二維空間以及不同流域的擴(kuò)展使廣義流態(tài)化相圖的繪制成為可能。根據(jù)Kwauk廣義流態(tài)化相圖的繪制思路79，遍歷不同操作范圍的氣體速度（Ug*=Ug/Ut）和顆粒速度(Ud*=Up/

36、Ut)，判斷其所處流域，進(jìn)而采用相應(yīng)的EMMS模型或擴(kuò)展模型計(jì)算床層的平均空隙率，最終可得如圖11所示的操作相圖。值得說(shuō)明的是，該相圖給出了典型流域轉(zhuǎn)變時(shí)的特征速度曲線(xiàn)，如噎塞曲線(xiàn)、液泛曲線(xiàn)和湍動(dòng)流態(tài)化速度的參考值。噎塞的預(yù)測(cè)作為EMMS模型的重要特征已經(jīng)被廣泛討論，這里不再贅述。而湍動(dòng)流域的確定一直存在廣泛爭(zhēng)議，不同研究者對(duì)湍動(dòng)氣速的確定標(biāo)準(zhǔn)持有不同的觀點(diǎn)。本研究根據(jù)基于團(tuán)聚物的原始EMMS模型和基于氣泡的鼓泡EMMS模型聯(lián)立求解，所獲得的二者的交叉點(diǎn)被認(rèn)為是連續(xù)相從密相轉(zhuǎn)變?yōu)橄∠嗟南噢D(zhuǎn)變點(diǎn)，并將其所對(duì)應(yīng)的流體速度定義為湍動(dòng)流態(tài)化速度，如圖中黑色虛線(xiàn)Uc所示。而液泛作為另外一個(gè)重要的流域轉(zhuǎn)折

37、特征，一般發(fā)生在當(dāng)逆流床的氣速或顆粒通量過(guò)大時(shí)，顆粒被氣流帶出從而中斷系統(tǒng)操作的狀態(tài)。通過(guò)考察EMMS逆流下行床模型解的特點(diǎn)，本研究將空隙率對(duì)氣速二階導(dǎo)的最大值點(diǎn)設(shè)為液泛發(fā)生速度，而EMMS無(wú)解點(diǎn)設(shè)為液泛終止速度65，并通過(guò)和實(shí)驗(yàn)對(duì)比驗(yàn)證了這種定義的合理性。圖11圖11Geldart B類(lèi)顆粒的EMMS廣義流態(tài)化相圖69Fig.11EMMS-based generalized fluidization phase diagrams for Geldart B particles693.2 循環(huán)流化床全回路穩(wěn)態(tài)建模及軟件開(kāi)發(fā)基于以上各擴(kuò)展模型和操作相圖，本課題組進(jìn)一步完善了基于EMMS方法的循環(huán)

38、流化床全回路穩(wěn)態(tài)計(jì)算模型70,72。圖12(a) 展示了循環(huán)流化床全回路建模的基本流程：（1）首先根據(jù)系統(tǒng)的邊界和接口情況將其劃分為不同的循環(huán)通路和相應(yīng)的子反應(yīng)器組合；（2）遍歷各循環(huán)通路內(nèi)的顆粒循環(huán)量等參數(shù)，并根據(jù)所處的氣固流速選擇合適的計(jì)算模型，如果操作處于廣義流態(tài)化狀態(tài)，則選取EMMS擴(kuò)展模型計(jì)算或通過(guò)相圖查閱該單元內(nèi)的動(dòng)力學(xué)參數(shù)，其他單元模塊如旋風(fēng)、出口管等則使用經(jīng)驗(yàn)關(guān)聯(lián)式求解；（3）根據(jù)系統(tǒng)的壓降平衡和物料平衡迭代計(jì)算直到收斂。圖12(b) 是一套多產(chǎn)異構(gòu)烷烴的催化裂化過(guò)程冷態(tài)實(shí)驗(yàn)裝置。圖12(c)是對(duì)該反應(yīng)器進(jìn)行全回路計(jì)算的結(jié)果，可見(jiàn)模型可以合理預(yù)測(cè)反應(yīng)器內(nèi)流動(dòng)參數(shù)的宏觀分布。圖1

39、2圖12基于EMMS的循環(huán)流化床統(tǒng)一建模72Fig.12EMMS-based general method for global hydrodynamic calculation of typical CFB systems72為實(shí)現(xiàn)流化床反應(yīng)器整體穩(wěn)態(tài)動(dòng)力學(xué)的快速計(jì)算，EMMS全回路建模方法還被封裝成軟件Virtual Fluidization（2022SRBJ0279）。軟件采用模塊化設(shè)計(jì)，可以根據(jù)實(shí)際的工業(yè)裝置自定義搭建、更改和優(yōu)化反應(yīng)器設(shè)計(jì)。軟件系統(tǒng)通過(guò)調(diào)用全循環(huán)計(jì)算代碼可以在幾秒到幾分鐘之內(nèi)完成穩(wěn)態(tài)動(dòng)力學(xué)的計(jì)算和顯示。如圖13所示，系統(tǒng)合理預(yù)測(cè)了不同操作條件下反應(yīng)器內(nèi)的各模塊的顆粒濃度非均勻分布。通過(guò)導(dǎo)出整體分布數(shù)據(jù)作為后續(xù)CFD計(jì)算的初場(chǎng)，可以提高CFD模擬的效率和速度。圖13圖13采用Virtual Fluidization軟件對(duì)循環(huán)流化床系統(tǒng)的穩(wěn)態(tài)動(dòng)力學(xué)