王利民， 付少童

(1.中國科學(xué)院過程工程研究所 多相復(fù)雜系統(tǒng)國家重點(diǎn)實(shí)驗(yàn)室，北京 100190；2.中國科學(xué)院大學(xué) 化學(xué)工程學(xué)院，北京 100049)

1 引 言

顆粒流體系統(tǒng)廣泛存在于自然界和工業(yè)過程，由于其涉及復(fù)雜的多相流動(dòng)過程和時(shí)空多尺度結(jié)構(gòu)，通過實(shí)驗(yàn)或理論分析存在較大困難，因此基于數(shù)值模擬展開顆粒流體系統(tǒng)的研究已成為不可或缺的重要手段。根據(jù)對流體相和顆粒相描述精度的不同，顆粒流體系統(tǒng)的模擬方法通常分為顆粒解析直接數(shù)值模擬PR-DNS(particle-resolved direct numerical simulation)、離散顆粒模擬DPS(discrete particle simulation)和雙流體模型TFM(two -fluid model)三類[1]。

顆粒流體系統(tǒng)模擬中,TFM的顆粒相處理為連續(xù)介質(zhì)的擬流體，PR-DNS和DPS的顆粒相通常是離散處理，根據(jù)牛頓第二定律對顆粒進(jìn)行 Lagrange 跟蹤。目前在顆粒流體系統(tǒng)的PR-DNS中，常忽略或簡單處理顆粒之間的作用。然而，大多數(shù)體系中顆粒間的相互作用至關(guān)重要，需要準(zhǔn)確和高效地對其進(jìn)行描述。顆粒碰撞處理主要有軟球模型、硬球模型和直接Monte Carlo方法DSMC(direct simulation Monte Carlo)。軟球模型的典型代表是離散單元法DEM(discrete element method)[2]，該方法基于彈性阻尼和滑移等機(jī)制計(jì)算顆粒間的碰撞力，再根據(jù)牛頓第二定律采用時(shí)間驅(qū)動(dòng)算法更新顆粒的速度和位置；硬球模型將顆粒碰撞處理成二體碰撞，其不求解顆粒受力，基于動(dòng)量守恒原理直接獲得顆粒碰撞后的速度信息；DSMC[3]同樣采用硬球模型更新碰撞后顆粒速度與角速度，但運(yùn)用概率抽樣確定顆粒是否碰撞，因此該方法也可視為硬球模型的一種變化形式。其中，軟球模型由于能夠獲得顆粒詳盡的受力信息而得到廣泛應(yīng)用。

顆粒流體系統(tǒng)中的流體相主要基于連續(xù)介質(zhì)假設(shè)框架下的N-S(Navier-Stokes)方程或其簡化形式，采用常規(guī)的有限差分、有限體積和有限元等方法進(jìn)行求解。由于該類方法流程復(fù)雜且包含隱式算法，在并行化和大規(guī)模計(jì)算上還存在困難[4]。格子Boltzmann方法LBM(lattice Boltzmann method)[5]是20世紀(jì)80年代末逐漸發(fā)展起來的一種流體系統(tǒng)建模和模擬的新方法。從離散的網(wǎng)格來說，其具有歐拉方法的屬性；從離散的粒子觀點(diǎn)來看，又具備拉格朗日方法的特點(diǎn)。這種特性使得該方法的研究思路與傳統(tǒng)的流體模擬方法完全不同，其屬于介于宏觀連續(xù)介質(zhì)模型和微觀分子動(dòng)力學(xué)模型之間的介觀模型，且物理背景清晰，受到國內(nèi)外學(xué)者的廣泛關(guān)注。隨著LBM基礎(chǔ)理論不斷完善，目前其研究已經(jīng)從理論向工程應(yīng)用發(fā)展，并在許多復(fù)雜物理現(xiàn)象和過程進(jìn)行了大量的應(yīng)用[6]，如湍流DNS、大渦模擬和耦合各種湍流模型的模擬等；針對多組分系統(tǒng)和不混溶的多流體系統(tǒng)，分別建立了顏色模型、偽勢模型、自由能模型和動(dòng)力學(xué)模型等；針對化學(xué)反應(yīng)系統(tǒng)，建立擴(kuò)散和燃燒模型。此外，LBM已拓展到模擬非牛頓流體、磁流體、聚合物流動(dòng)、血液流、光子運(yùn)動(dòng)、聲場和微流動(dòng)等復(fù)雜體系[6]。

工業(yè)流化床反應(yīng)器的最大流動(dòng)結(jié)構(gòu)可達(dá)到米量級，但這些結(jié)構(gòu)卻受發(fā)生在毫米量級的顆粒-顆粒間碰撞和顆粒-流體間相互作用的直接影響。目前很難通過某一具體的模擬方法涵蓋所有的時(shí)間和空間尺度，根據(jù)關(guān)注重點(diǎn)問題的不同，應(yīng)用格子Boltzmann框架可在不同時(shí)空尺度和精度上開展顆粒流體系統(tǒng)模擬(圖1)，如在微觀層次，LB -based PR-DNS擬獲得的詳細(xì)信息深入探索顆粒流體界面的流動(dòng)、傳遞和反應(yīng)以及兩相相互作用的本構(gòu)關(guān)系；在介觀層次，利用 LB -based DPS探索實(shí)驗(yàn)室規(guī)模的顆粒流體流動(dòng)，建立更為準(zhǔn)確的宏觀多相流動(dòng)模型；在宏觀層次，發(fā)展格子Boltzmann雙流體模型(LB -based TFM)有望實(shí)現(xiàn)工業(yè)級別的運(yùn)算。

圖1 基于格子Boltzmann的顆粒流體系統(tǒng)計(jì)算框架

本文將重點(diǎn)介紹LBM在顆粒流體系統(tǒng)模擬方面的研究進(jìn)展，包括 LB -based PR-DNS方法、LB -based DPS方法及其最新應(yīng)用，以及關(guān)于開發(fā) LB -based TFM研究的一些進(jìn)展。

2 顆粒解析直接數(shù)值模擬

顆粒流體系統(tǒng)PR-DNS與單相湍流DNS相類似，其中流體相的控制方程直接采用數(shù)值計(jì)算求解，無任何湍流模型，流場中顆粒周圍計(jì)算網(wǎng)格縮小到顆粒尺度以下，能夠分辨包含Kolmogorov耗散尺度以內(nèi)的所有空間尺度，而顆粒的受力則通過對其表面的黏性力與壓力進(jìn)行積分獲得，不引入經(jīng)驗(yàn)?zāi)Ｐ蚚7]。根據(jù)流體相的描述方式不同，顆粒流體系統(tǒng)的PR-DNS又可以分成三類(圖2)，一是基于N-S方程的PR-DNS方法；二是基于格子的PR-DNS方法；三是基于粒子的PR-DNS方法。LB -based PR-DNS方法屬于基于格子的PR-DNS方法。根據(jù)流固作用處理方式不同，基于LBM的顆粒流體系統(tǒng)PR-DNS分為邊界鏈法和格點(diǎn)法。

圖2 顆粒流體系統(tǒng)的顆粒解析直接數(shù)值模擬

2.1 邊界鏈法

Ladd[8]率先提出了基于有限體積顆粒的LBM模擬液固懸浮，通過修正的回彈邊界保證運(yùn)動(dòng)邊界的無滑移，通過動(dòng)量交換計(jì)算流體和固體顆粒之間的作用力。Ladd方法中固體顆粒實(shí)際是由一個(gè)階梯狀的殼組成，固體顆粒邊界像許多流體-固體格子鏈組成的，因此稱為邊界鏈法。

Ladd方法[8]的流固作用通過流體粒子和固體顆粒之間的動(dòng)量交換實(shí)現(xiàn)，其優(yōu)點(diǎn)是計(jì)算直觀簡潔，而且保證了流體的質(zhì)量和動(dòng)量局部守恒。但這種階梯狀的殼近似表示真實(shí)固體顆粒，往往會(huì)引入誤差，而且在大多數(shù)情況下，這種誤差會(huì)導(dǎo)致顆粒的動(dòng)力學(xué)半徑略大于其真實(shí)半徑。只有在固體顆粒半徑遠(yuǎn)大于LBM網(wǎng)格時(shí)，固體顆粒的動(dòng)力學(xué)半徑才和其真實(shí)半徑相近。此外，Ladd方法中，固體顆粒內(nèi)部由流體粒子占據(jù)，要求固體與流體密度比大于1。針對這一情況，Aidun等[9]提出不使用顆粒內(nèi)部流體的方法，能夠處理固體密度小于流體密度的問題。在Ladd方法的基礎(chǔ)上，有學(xué)者進(jìn)一步提出了滿足局部質(zhì)量守恒的方法[10]，并模擬了多顆粒沉降問題。然而，以Ladd為基礎(chǔ)的方法均存在著邊界形狀與真實(shí)邊界并不相符的問題，為了解決階梯狀近似，不少學(xué)者提出了改進(jìn)方法，大部分都是采用處理曲面邊界的方法,如曲面邊界插值[11]。而這些處理方法雖然提高了邊界處理的精確性，但在計(jì)算效率上卻付出了巨大的代價(jià)，通常處理顆粒數(shù)目較少。但基于邊界鏈法已可實(shí)現(xiàn)模擬含17845個(gè)顆粒的湍管流[12]。

2.2 格點(diǎn)法

Feng等[13]拋棄了固定顆粒邊界的想法，把顆粒的邊界當(dāng)作可變形，結(jié)合浸入邊界法，提供了一種新的處理運(yùn)動(dòng)流固邊界的方法，這種方法避開了顆粒邊界無法精確描述的困難，但對顆粒表面剛性系數(shù)的選擇具有隨意性。Noble等[14]提出了浸入運(yùn)動(dòng)邊界法來處理流固耦合，通過在LBM演化方程中加入附加碰撞項(xiàng)，實(shí)現(xiàn)了對固體顆粒邊界相對準(zhǔn)確的描述，同時(shí)沿用了LBM的計(jì)算體系，保留了原算法的優(yōu)點(diǎn)。Noble等在浸入運(yùn)動(dòng)邊界法中引入了格子控制體和格子固含率的概念。格子固含率為格子控制體內(nèi)顆粒體積占格子控制體體積的比重，固含率越高表示格子內(nèi)包含的顆粒組分越多，相應(yīng)的流體格子與顆粒的相互作用越強(qiáng)。將格子固含率作為重要參數(shù)進(jìn)一步得出權(quán)重系數(shù),實(shí)現(xiàn)對顆粒邊界的光滑描述。權(quán)重系數(shù)存在多種構(gòu)建形式[15-17]，均在一定程度上改進(jìn)了模擬的精度。

Cook等[18]耦合浸入運(yùn)動(dòng)邊界法和DEM模擬了顆粒相互作用占主導(dǎo)的顆粒流。馮云田等[19]在此基礎(chǔ)上，引入大渦模擬LES(large eddy simulation)，成功地對高雷諾數(shù)顆粒流體系統(tǒng)進(jìn)行了模擬。王利民等[20]發(fā)展了一種基于粒子-格子耦合方法的PR-DNS方法，該方法中顆粒模型采用結(jié)合硬球模型和軟球模型各自優(yōu)點(diǎn)的時(shí)間驅(qū)動(dòng)硬球模型TDHS(time -driven hard-sphere)[21]進(jìn)行求解，該算法將顆粒/流體密度比從2.5提高到了1500，成功復(fù)現(xiàn)了散式流態(tài)化和聚式流態(tài)化的典型現(xiàn)象[20]。在此基礎(chǔ)上，周國峰等[22]成功模擬了二維和三維DKT(Drafting-Kissing-Tumbling)過程，驗(yàn)證了該算法的有效性，此外，還實(shí)現(xiàn)了Lee-Edwards邊界條件與基于格點(diǎn) LB -based PR-DNS的耦合[23]，克服了邊界鏈接法無法滿足伽利略不變性的問題，并成功預(yù)測了低雷諾數(shù)下顆粒懸浮體系的表觀粘度。熊勤鋼等[7]借鑒基于粒子-格子耦合方法的顆粒流體系統(tǒng)PR-DNS的思想，使用簡化DEM代替硬球算法，利用中國科學(xué)院過程工程研究所研制的千萬億次超級計(jì)算系統(tǒng) Mole -8.5，使用672個(gè)GPU進(jìn)行了大規(guī)模并行計(jì)算，實(shí)現(xiàn)了全球最大規(guī)模顆粒流體系統(tǒng)PR-DNS，其中二維模擬包含1166400個(gè)顆粒，三維模擬包含129024個(gè)顆粒(圖3)。

圖3 顆粒懸浮系統(tǒng)大規(guī)模直接模擬粒子和流場分布[7]

周國峰等[24]分析了非均勻氣固系統(tǒng)中顆粒受到的曳力情況，對氣固流態(tài)化系統(tǒng)進(jìn)行了三維大規(guī)模PR-DNS，成功復(fù)現(xiàn)了流動(dòng)的非均勻結(jié)構(gòu)，揭示了沒有考慮非均勻結(jié)構(gòu)的Wen &Yu曳力公式預(yù)測偏高，且受力方向存在偏差。

劉曉雯等[25,26]在中國科學(xué)院過程工程研究所研制的千萬億次超級計(jì)算系統(tǒng) Mole -8.5E 上，使用192個(gè)Tesla K80 GPU運(yùn)行17天，對4.83×109個(gè)流體網(wǎng)格和115200個(gè)固體顆粒的周期懸浮進(jìn)行了直接數(shù)值模擬(圖4)，探索了介尺度結(jié)構(gòu)對顆粒統(tǒng)計(jì)屬性、氣固相間作用力以及傳遞特性的影響。

圖4 大規(guī)模氣固懸浮系統(tǒng)直接模擬瞬態(tài)結(jié)果[26]

劉曉雯等[25,26]發(fā)現(xiàn)顆粒脈動(dòng)速度和局部平均無量綱曳力具備很強(qiáng)的尺度依賴性，提出了一種雙峰分布形式以修正顆粒脈動(dòng)速度分布，闡明了顆粒脈動(dòng)速度偏離Maxwell分布的內(nèi)在機(jī)理，并從微觀層次上揭示了介尺度曳力的尺度依賴性，提出了以Fr數(shù)反映結(jié)構(gòu)效應(yīng)的曳力關(guān)聯(lián)式[26]。

此外，耦合直接力浸入邊界法的LBM應(yīng)用于三維顆粒懸浮PR-DNS，如1024個(gè)熱顆粒在流化床中的運(yùn)動(dòng)[27]；Rubinstein等[28]模擬了2000個(gè)顆粒在流化床中的運(yùn)動(dòng)；Eshginejadfard等[29]對305個(gè)顆粒槽道流進(jìn)行了模擬。LB -based PR-DNS也應(yīng)用于分析多顆粒沉降[30,31]或顆粒自旋及振動(dòng)對流動(dòng)的影響[32-34]、圓柱與柔性體耦合[35]、管道泄漏時(shí)土壤的流化現(xiàn)象[36]、泥石流[37]以及耦合傳熱的顆粒流問題[38]。通過對LBM中的弛豫時(shí)間進(jìn)行修正，LB -based PR-DNS擴(kuò)展到研究非牛頓流體顆粒流體系統(tǒng)中的顆粒流動(dòng)規(guī)律[39-44]。通過重構(gòu)顆粒的形式和顆粒間作用力，基于浸入邊界法的LBM可以非常容易地處理諸如多邊形顆粒[45]、可變形顆粒[46]和粘性顆粒[47]的遷移等，甚至可以分析巖土工程中的水力壓裂問題[48]。以上均顯示了基于 LB -based PR-DNS的計(jì)算策略在模擬顆粒流體系統(tǒng)中的巨大優(yōu)勢和潛力。

3 離散顆粒模擬

離散顆粒模擬(或顆粒軌道模型)將流體相處理為連續(xù)介質(zhì)，固相處理為離散顆粒。在歐拉坐標(biāo)系考察流體運(yùn)動(dòng)，在拉格朗日坐標(biāo)系以及顆粒層次分析固相運(yùn)動(dòng)[49]。由流固耦合的實(shí)現(xiàn)方式不同，離散顆粒模擬一般可分三類，(1)單向耦合(one-way coupling)，只考慮流體對顆粒的作用，不考慮顆粒對流體的影響；(2)雙向耦合(two-way coupling)，考慮流體和顆粒之間的相互作用；(3)四向耦合(four-way coupling)，除了考慮流體和顆粒之間的相互作用外，還考慮顆粒間碰撞?，F(xiàn)在通常指的離散顆粒模擬，都是四向耦合的確定性顆粒軌道模型，在考慮顆粒間相互作用的情況下，采用曳力公式對顆粒與流體的相互作用進(jìn)行描述。同時(shí)，基于牛頓第三定律可實(shí)現(xiàn)流場與顆粒的相互關(guān)聯(lián)，即在各個(gè)控制體中，顆粒受到流體的作用力等于顆粒對流體施加的作用力。

在傳統(tǒng)的離散顆粒模擬中，流體相的控制方程是基于體積平均Navier-Stokes (volume-averaged Navier-Stokes,VANS)方程或者Navier-Stokes方程的簡化形式，通常采用SIMPLE(semi-implicit method for pressure-linked equations)算法等求解，這些處理方式都不易實(shí)現(xiàn)并行化,難以應(yīng)用于大規(guī)模的計(jì)算中[4]。基于LBM的算法具有良好的并行性，王利民等[1,50]將LBM應(yīng)用于離散顆粒模擬中流體相的流動(dòng)(圖5)，提出了基于格子 Boltzmann 的離散顆粒模擬( LB -based DPS)，其使用一種既考慮孔隙率又考慮氣固相對滑移速度的修正LBM求解流體流動(dòng)，取代了基于VANS方程的求解傳統(tǒng)，固體運(yùn)動(dòng)采用TDHS求解，即固體顆粒的碰撞為硬球碰撞，在固定的時(shí)間步長內(nèi)檢索顆粒間是否發(fā)生碰撞以及更新顆粒運(yùn)動(dòng)。由于在顆粒濃度較高或顆粒溫度較低的場景下應(yīng)用硬球模型還有一定的應(yīng)用局限性，王利民等[51]進(jìn)一步發(fā)展了一種介尺度LBM-DEM耦合模型。熊勤鋼等[52]也提出了一種快速的LBM-DEM耦合模型，但對力源項(xiàng)的添加方式有所不同。

圖5 D2Q9計(jì)算流體格點(diǎn)和固體顆粒[1]

(1)

式中v為格點(diǎn)處的流體速度，vs為顆粒在當(dāng)前格點(diǎn)的平均速度，-i表示與i方向相反。修正后的格子Boltzmann方程為[1]

(2)

式中εs和εg分別為格子控制體內(nèi)部固體顆粒和氣體的體積分率，且εs+εg=1。

LB -based DPS成功模擬了二維單孔射流鼓泡床[1,50,52]和三維鼓泡床(圖6)的流態(tài)化過程，得到了與經(jīng)驗(yàn)關(guān)聯(lián)式和實(shí)驗(yàn)數(shù)據(jù)相吻合的模擬結(jié)果。李斌等[53,54]將 LB -based DPS進(jìn)一步應(yīng)用于二維噴動(dòng)床和鼓泡床的流化過程，獲得了與實(shí)驗(yàn)數(shù)據(jù)相吻合的結(jié)果。以上均驗(yàn)證了 LB -based DPS能夠得到流固系統(tǒng)更多的細(xì)節(jié)信息和內(nèi)在特性，此外，LB -based DPS在計(jì)算精度、耗費(fèi)時(shí)間和效率之間能達(dá)到很好的平衡,可以直接擴(kuò)展到與其他固相求解器耦合計(jì)算。

圖6 三維流化床離散顆粒模擬流動(dòng)結(jié)構(gòu)的演化

4 雙流體模型

目前 LB -based TFM研究文獻(xiàn)較少，傳統(tǒng)LBM的顏色模型、偽勢模型、自由能模型和動(dòng)力學(xué)模型只適用于不混溶的二元流體模擬，不能應(yīng)用于流體相互滲透的多相流雙流體層次模擬。王鐵峰等[55]率先開展利用LBM構(gòu)建雙流體模型的研究，其利用無量綱方法變換格子Boltzmann方程并添加壓力修正項(xiàng)來實(shí)現(xiàn)，但其恢復(fù)的宏觀方程中，密度為兩相的合密度而非分相密度；Sankaranarayanan等[56]從連續(xù)Boltzmann方程出發(fā)，考慮顆粒溫度方程，提出了一種隱式格子Boltzmann模型，此模型也沒有真正解決問題。實(shí)現(xiàn) LB -based TFM的關(guān)鍵在于構(gòu)造能恢復(fù)為VANS方程的格子Boltzmann模型，郭照立等[57]構(gòu)造了一個(gè)可用于求解VANS方程的格子Boltzmann模型，該模型的平衡分布函數(shù)和作用力項(xiàng)都包含孔隙率，作用力的計(jì)算需要給定組成多孔介質(zhì)的顆粒粒徑，適合模擬均勻孔隙率的多孔介質(zhì)內(nèi)部流動(dòng)，對于具有孔隙率梯度的復(fù)雜多相流問題，顯式添加的離散力項(xiàng)通常導(dǎo)致模型不易收斂。Sungkon等[58]將VANS方程中帶有體積分?jǐn)?shù)的項(xiàng)均移到方程右端作為力源項(xiàng)，變換得到帶有源項(xiàng)的N-S方程，并構(gòu)造能恢復(fù)帶有源項(xiàng)的N-S方程的LB方程，但是所得LB方程只有一階精度。宋飛飛等[59]類似地將VANS方程轉(zhuǎn)為帶有源項(xiàng)的N-S方程，利用帶有源項(xiàng)的作用力LBM模型發(fā)展了具有二階精度的隱式格子Boltzmann方程，處理較復(fù)雜，不易并行。王利民等[1]考慮了體積分?jǐn)?shù)和相間滑移速度對流體相流動(dòng)的影響，基于浸入邊界的思想，提出了修正的格子Boltzmann方程，用于描述氣相流動(dòng)，嘗試拓展到雙流體模擬。張金鳳等[4]提出基于VANS方程的格子Boltzmann模型，通過Chapman-Enskog分析證明提出模型確實(shí)能夠恢復(fù)到VANS方程。Blais等[60]指出體系存在孔隙率梯度的情況下，直接添加力源項(xiàng)修正壓力梯度的LBM形式是不穩(wěn)定的，并通過重構(gòu)LBM中的碰撞算子以還原成VANS方程，其在大孔隙率梯度的情況下具備更好的穩(wěn)定性。H?cker等[61]則是通過修正遷移步驟以提高穩(wěn)定性，結(jié)合對流擴(kuò)散的LBM實(shí)現(xiàn)了Rayleigh-Taylor不穩(wěn)定性的模擬，但如何解決TFM下顆粒的緊密堆積還需進(jìn)一步研究。雖然存在多種思路構(gòu)建VANS方程的格子Boltzmann模型，但在使用上均存在一定局限性，將其應(yīng)用到兩相流的研究仍需進(jìn)一步探索。

5 格子多相流體力學(xué)模擬軟件

LB -based PR-DNS和DPS基于的浸入運(yùn)動(dòng)邊界法以固含率將流體、顆粒和邊界統(tǒng)一描述，保留了LB求解的并行性和遷移的簡單性；同時(shí)格子Boltzmann方程基于笛卡爾網(wǎng)格求解能夠處理復(fù)雜邊界，避免了傳統(tǒng)CFD繪制復(fù)雜構(gòu)體網(wǎng)格的困難，在大規(guī)模并行條件下仿真規(guī)模已具備由實(shí)驗(yàn)室尺度向工業(yè)尺度的可行性?；诖?，中國科學(xué)院過程工程研究所EMMS團(tuán)隊(duì)集成相關(guān)算法開發(fā)了以顆粒流仿真和多GPU并行為特色的高性能格子多相流體力學(xué)LMFD(lattice-based multi-fluids dynamics)模擬軟件，如圖7所示。LMFD基于C/C++開發(fā)，目前穩(wěn)定版本支持二維 LB -based PR-DNS和DPS，以及三維簡單構(gòu)型的 LB -based PR-DNS計(jì)算。

圖7 LMFD軟件界面

顆粒流體仿真的主流商業(yè)軟件有Barracuda，F(xiàn)luent和開源的軟件如MFiX和OpenFOAM。以上產(chǎn)品均基于傳統(tǒng)的有限體積法開發(fā)，在并行性上受到很多制約。而LB基的主流求解器，如PowerFlow，XFlow和Palabos均以求解單相流為主，而以多相流仿真為主的LBM求解器還比較少。Seil等[62]耦合了Palabos與LIGGGHTS兩個(gè)開源軟件，實(shí)現(xiàn)了 LB -based PR-DNS求解器LBDEMcoupling。G?tz等[63]在具備自適應(yīng)網(wǎng)格加密的單相流程序WaLBerla基礎(chǔ)上實(shí)現(xiàn)了 LB -based PR-DNS，用于研究河床中的沙丘形成。但以上程序的應(yīng)用場景較為有限，且均無法實(shí)現(xiàn) LB -based DPS。LMFD軟件通過構(gòu)建通用的LB基高性能氣固兩相流離散模擬平臺(tái)，再逐步形成完備的多相流求解模塊，將有望彌補(bǔ)這一空缺，在化工、冶金、能源、水利、環(huán)保和航空航天等相關(guān)領(lǐng)域具有廣泛的應(yīng)用前景。

6 結(jié) 論

LBM作為一種近幾十年發(fā)展的流體系統(tǒng)建模及模擬新方法，為顆粒流體系統(tǒng)的建模及計(jì)算提供了一種新思路。本文對LBM應(yīng)用于顆粒流體系統(tǒng)的進(jìn)展作了較為系統(tǒng)總結(jié)，針對不同層次的顆粒流體計(jì)算框架，未來重點(diǎn)的研究方向如下。

(1) 對于格子Boltzmann顆粒解析直接數(shù)值模擬，隨著計(jì)算機(jī)硬件的發(fā)展，高雷諾數(shù)顆粒懸浮、顆粒湍流作用、復(fù)雜顆粒、耦合傳熱、傳質(zhì)及化學(xué)反應(yīng)等復(fù)雜過程等基礎(chǔ)研究將會(huì)受到廣泛關(guān)注。

(2) 對于格子Boltzmann離散顆粒模擬，擴(kuò)大計(jì)算規(guī)模應(yīng)用于工業(yè)反應(yīng)器是未來發(fā)展趨勢，如顆粒相采用多相流網(wǎng)格質(zhì)點(diǎn)群(MP-PIC)方法中parcel的處理方式，或采用GPU等硬件加速計(jì)算等。

(3) 對于格子Boltzmann雙流體模擬，隨著LBM基礎(chǔ)理論研究深入開展，LBM代替?zhèn)鹘y(tǒng)有限元或有限體積方法求解雙流體模型方程在不久的將來肯定會(huì)得以實(shí)現(xiàn)。

致謝:感謝李靜海院士對顆粒多相流模擬工作的長期支持與指導(dǎo)；感謝葛蔚研究員和王小偉研究員在格子Boltzmann模擬方面長期友好合作；感謝曾參與本文涉及工作的熊勤鋼、周國峰、劉曉雯、魏民、張博和張金鳳等畢業(yè)同學(xué)；感謝陳建華博士對本文工作的積極建議和幫助。