張瀟丹，雍玉梅，李文軍，趙元生，李媛媛，楊巧文，楊超

（1中國(guó)礦業(yè)大學(xué)（北京）化學(xué)與環(huán)境工程學(xué)院，北京 100083；2中國(guó)科學(xué)院過程工程研究所綠色過程與工程重點(diǎn)實(shí)驗(yàn)室，北京 100190；3華北科技學(xué)院環(huán)境工程學(xué)院，河北 廊坊 065201；4中國(guó)石油化工研究院渣油加氫實(shí)驗(yàn)室，北京 102200）

?

REV尺度多孔介質(zhì)格子Boltzmann方法的數(shù)學(xué)模型及應(yīng)用的研究進(jìn)展

張瀟丹1，2，雍玉梅2，李文軍3，趙元生4，李媛媛2，楊巧文1，楊超2

（1中國(guó)礦業(yè)大學(xué)（北京）化學(xué)與環(huán)境工程學(xué)院，北京 100083；2中國(guó)科學(xué)院過程工程研究所綠色過程與工程重點(diǎn)實(shí)驗(yàn)室，北京 100190；3華北科技學(xué)院環(huán)境工程學(xué)院，河北 廊坊 065201；4中國(guó)石油化工研究院渣油加氫實(shí)驗(yàn)室，北京 102200）

摘要：綜述了多孔介質(zhì)表征體元尺度（REV）格子Boltzmann模型的研究進(jìn)展，根據(jù)對(duì)多孔介質(zhì)處理方式主要分為部分反彈模型和阻力模型兩類，分析歸納了各類模型的優(yōu)缺點(diǎn)。由于阻力模型中滲流的廣義格子Boltzmann方程（GLBE）的作用力是基于GUO等的作用力模型，可以準(zhǔn)確得到宏觀方程，不存在離散誤差，且模型的平衡分布函數(shù)和作用力項(xiàng)中都包含反應(yīng)介質(zhì)特性的孔隙率，因而應(yīng)用最為廣泛。本文還重點(diǎn)介紹了REV尺度多孔介質(zhì)LBE模型在流動(dòng)、傳熱、傳質(zhì)、化學(xué)反應(yīng)及相變等過程中的具體應(yīng)用，認(rèn)為REV尺度多孔介質(zhì)內(nèi)的三傳一反數(shù)學(xué)模型中需要加入孔隙尺度因素，在更大工程尺度上應(yīng)該考慮過程參數(shù)的各向異性，展望了REV尺度多孔介質(zhì)LBE模型的發(fā)展和應(yīng)用前景。

關(guān)鍵詞：多孔介質(zhì)；表征體元尺度；格子Boltzmann方法；流動(dòng)；傳熱；傳質(zhì)

第一作者：張瀟丹（1989—），女，碩士研究生，主要從事化學(xué)工程數(shù)值模擬。聯(lián)系人：雍玉梅，研究員，主要從事化學(xué)工程數(shù)值模擬E-mail ymyong@ipe.ac.cn。

流體在多孔介質(zhì)內(nèi)的流動(dòng)是許多領(lǐng)域都涉及的一類復(fù)雜流動(dòng)，在能源化工、生物醫(yī)藥、空氣凈化、污水處理等領(lǐng)域都有廣泛應(yīng)用。多孔介質(zhì)內(nèi)的流動(dòng)是一種典型的多尺度問題，通常涉及3個(gè)尺度，即孔隙尺度、表征體元尺度（REV）和宏觀尺度。格子Boltzmann方法（LBM）是介于宏觀和微觀兩種方法之間的介觀數(shù)值模擬方法，它既克服了宏觀方法需要處理復(fù)雜邊界的劣勢(shì)，又克服了微觀方法計(jì)算尺度受限的局限性，在模擬多孔介質(zhì)方面具有獨(dú)特優(yōu)勢(shì)[1]。由于LBM算法簡(jiǎn)單，計(jì)算網(wǎng)格容易生成，可以比較方便處理流體與邊界之間、不同流體組分之間等復(fù)雜的相互作用，且程序代碼簡(jiǎn)單容易并行，應(yīng)用逐漸廣泛。格子Boltzmann方法模擬多孔介質(zhì)滲流主要有孔隙尺度和REV尺度?？紫冻叨萀BM的研究對(duì)象是單個(gè)或者若干個(gè)孔隙內(nèi)的流體，把介質(zhì)的骨架作為流場(chǎng)的邊界，用來探索滲流機(jī)理和基本規(guī)律。REV尺度LBM是以一個(gè)控制體為研究對(duì)象，求解REV體積平均的宏觀流動(dòng)量，流體與介質(zhì)的相互作用通過適當(dāng)?shù)奶幚矸椒▉砻枋?。REV尺度LBM方法相對(duì)于孔隙尺度LBM方法，優(yōu)點(diǎn)是前者不需要介質(zhì)的細(xì)節(jié)結(jié)構(gòu)、僅僅依賴于介質(zhì)的統(tǒng)計(jì)參數(shù)，因而計(jì)算效率較高，比較適合大區(qū)域的工程滲流計(jì)算。已有文獻(xiàn)總結(jié)了LBM方法的基本要素[2]，介紹了LBM的起源、連續(xù)Boltzmann方程、邊界條件和數(shù)值穩(wěn)定性，討論了LBM在流體流動(dòng)、多相流和多組分流、反應(yīng)擴(kuò)散中的應(yīng)用。CHEN等[3]則對(duì)LBM在湍流、多孔介質(zhì)流、燃燒、多孔介質(zhì)內(nèi)燃燒的應(yīng)用進(jìn)行了綜述，李楠等[4]分析了LBM在模擬這些復(fù)雜流體問題時(shí)存在的優(yōu)缺點(diǎn)。SUKOP 等[5]通過比較兩種REV尺度LBE模型，提出Dardis 和McCloskey的LBE模型不能用于各向異性的多孔介質(zhì)的觀點(diǎn)，但沒有詳細(xì)論述兩模型的具體差異。目前還沒有REV尺度多孔介質(zhì)LBM方法的綜述，因而本文總結(jié)了REV尺度多孔介質(zhì)LBE各類數(shù)學(xué)模型，分析了各模型的優(yōu)缺點(diǎn)，綜述了該方法在REV尺度多孔介質(zhì)內(nèi)流動(dòng)、傳熱、傳質(zhì)、化學(xué)反應(yīng)和相變等領(lǐng)域的應(yīng)用，分析了未來應(yīng)用前景，并對(duì)其發(fā)展趨勢(shì)提出建議。

1　REV尺度LBE模型

REV尺度LBE模型根據(jù)對(duì)多孔介質(zhì)的處理方式不同主要分為兩類：①部分反彈模型[6-11]，此種模型對(duì)經(jīng)典的固壁回彈模型進(jìn)行修改來擬合達(dá)西公式；②阻力模型[12-21]，此模型通過施加外力公式，得到Stokes或Brinkman方程。

1.1 部分反彈模型

部分反彈模型由3個(gè)步構(gòu)成，如式(1)。前兩個(gè)步是碰撞步和流動(dòng)步，與標(biāo)準(zhǔn)LBE模型相同。第三個(gè)步是部分反彈步，該過程中流體分布函數(shù)根據(jù)該點(diǎn)的固體密度重新定義。部分反彈模型中引入了多孔介質(zhì)連續(xù)表示方法，并將反彈邊界的作用耦合到演化方程中[7]。

式中，i是速度離散方向；t是離散時(shí)間；f為分布函數(shù)；τ為松弛時(shí)間；Δi(x,t)是流體粒子與固體骨架之間的碰撞，表示多孔介質(zhì)對(duì)流體影響。

最早的部分反彈模型是GAO和SHARMA[6]提出的多孔介質(zhì)格子氣自動(dòng)模型，可延伸到LBM中（GS-LBE模型)。該模型能夠模擬非均質(zhì)多孔介質(zhì)內(nèi)的流動(dòng)。該模型中?i(x,t)如式(2)。

式中，ns是固體格點(diǎn)的平均散射密度，取值范圍為0～1，代表多孔介質(zhì)的固體結(jié)構(gòu)；是i的相反方向。

但是GS-LBE模型存在格子氣的統(tǒng)計(jì)噪聲。1998年，DARDIS和MCCLOSKEY[7]則提出了第一個(gè)部分反彈LBE模型（DM-LBE模型)，是對(duì)Brinkman方程的間接求解。與GS-LBE模型相比，DM-LBE模型用ns代替孔隙尺度多孔介質(zhì)的結(jié)構(gòu)參數(shù)布爾變量，消除了格子氣的不穩(wěn)定性。滲透率依賴于固體格點(diǎn)的平均散射密度ns，并且對(duì)應(yīng)的有效黏性系數(shù)與流體的黏性系數(shù)相同。該模型中為?i(x,t)如式(3)。

THORNE和SUKOP[8]根據(jù)部分反彈形式不同提出了第二類部分反彈格式（TS-LBE模型），該模型流出分布函數(shù)反映臨近點(diǎn)流入分布函數(shù)碰撞過程，該類反彈格式需計(jì)算每個(gè)格點(diǎn)的碰撞步來估算流出分布函數(shù)密度。該模型中?i(x,t)如式(4)。

式中t*是碰撞前的時(shí)間，滿足t＜t*＜t+1，fic為碰撞后的分布函數(shù)。

DM-LBE模型和TS-LBE模型的部分反彈格式需要根據(jù)臨近格點(diǎn)計(jì)算得到，計(jì)算復(fù)雜，不利于并行計(jì)算。WALSH、BURWINKLE和SAAR[9]引入了新的部分反彈方法（WBS-LBE模型），該方法在碰撞步前改變流入流體方向，同時(shí)計(jì)算碰撞步和反彈步，且計(jì)算中不包含臨近格點(diǎn)的流體分布函數(shù)，使計(jì)算更容易并行，從而節(jié)省計(jì)算時(shí)間。與DM-LBE模型和TS-LBE模型的部分反彈方法相比較，該方法在非均質(zhì)多孔介質(zhì)流動(dòng)過程中滿足質(zhì)量守恒，且提高了LB方法數(shù)值模擬浮力驅(qū)動(dòng)及應(yīng)力和濃度梯度引發(fā)的擴(kuò)散過程的精度。該模型中?i(x,t)如式(5)。

WBS-LBE模型雖然在非均質(zhì)多孔介質(zhì)內(nèi)滿足質(zhì)量守恒，更加適用于浮力驅(qū)動(dòng)的流動(dòng)和擴(kuò)散過程，但是該模型在不同滲透特性的多孔介質(zhì)界面處速度不連續(xù)，ns接近1時(shí)會(huì)導(dǎo)致數(shù)值不穩(wěn)定。YEHYA、NAJI和SUKOP[10]修正了該模型，提出新的模型（YNS-LBE模型）。該模型也在碰撞前改變流入流體方向，但是公式中使用i方向的分布函數(shù)代替i方向的分布函數(shù)，使其在界面處速度連續(xù)，如式(6)。

WBS-LBE模型在i方向上的流體分布函數(shù)是碰撞后的分布函數(shù)，而方向上的分布函數(shù)則是在碰撞前的分布函數(shù)，導(dǎo)致反彈回的流體不僅沿其反方向流體，也會(huì)向其他方向流動(dòng)，引起數(shù)值不穩(wěn)定性。ZHU和MA[11]重新分配同一時(shí)刻同一格點(diǎn)碰撞后方向相反的一對(duì)流體質(zhì)點(diǎn)，每個(gè)方向上碰撞后的流體質(zhì)點(diǎn)分為兩部分，一部分沿該方向流向下一格點(diǎn)，另一部分由于固體壁面的作用沿相反方向反彈回去，提出新的部分反彈模型（ZM-LBE模型），該模型在同一時(shí)間的格點(diǎn)僅改變方向相反的一對(duì)流體質(zhì)點(diǎn)，因此該模型滿足質(zhì)量守恒，且數(shù)值穩(wěn)定性比WBS-LBE模型好。該模型的松弛時(shí)間接近1時(shí)，有效黏度系數(shù)與流體黏度系數(shù)相等。該模型中?i(x,t)如式(7)。

部分反彈模型區(qū)別主要在于如何處理反彈邊界的作用、有效黏度、平均散射密度及有效滲透率，總結(jié)目前主要部分反彈模型見表1。表中包括DM-LBE模型、TS-LBE模型、WBS-LBE模型、YNS-LBE模型和ZM-LBE模型，這5個(gè)模型的滲透率根據(jù)每個(gè)格點(diǎn)的平均散射密度ns計(jì)算得到。前4個(gè)模型在ns=0時(shí)，均可簡(jiǎn)化為標(biāo)準(zhǔn)LBE模型. 在ns=0時(shí)，ZM-LBE模型不能自動(dòng)滿足反彈邊界，因?yàn)樵撃Ｐ驮谕粫r(shí)間格點(diǎn)僅改變方向相反的一對(duì)流體質(zhì)點(diǎn)。WBS-LBE模型、YNS-LBE模型和ZM-LBE模型在ns=1時(shí)能夠滿足全反彈LBE模型，其余模型不滿足全反彈LBE模型，因?yàn)楫?dāng)ns=1時(shí)，DM-LBE部分反彈模型在不透水邊界處有質(zhì)量損失，而TS-LBE部分反彈格式模型僅反映一個(gè)時(shí)間步的流體函數(shù)而不是原則上的兩個(gè)時(shí)間步。其次，GS-LBE模型、WBS-LBE模型、YNS-LBE模型和ZM-LBE模型能夠滿足質(zhì)量守恒，DM-LBE模型和TS-LBE模型在高滲透性和低滲透性多孔介質(zhì)邊界之間不對(duì)稱而不滿足質(zhì)量守恒，需要在格點(diǎn)的邊緣設(shè)置ns來修正模型，只有在密度恒定或者平均散射密度不變時(shí)才滿足質(zhì)量守恒。CHEN等[22]比較了GS-LBE模型、DM-LBE模型和TS-LBE模型，分析了在松弛時(shí)間τ=1時(shí)，為了保證流體黏度為正值，3個(gè)模型中ns的上限值。GS-LBE模型中ns< (1/2)≈0.289 DM-LBE模型中ns< (–1/28)≈0.234，TS-LBE模型和WBS-LBE模型中ns在0～1之間。ZM-LBE模型的ns上限為0.5，說明ZM-LBE模型最多有一半的流體質(zhì)點(diǎn)可以反彈回相反的方向。綜上，YNS-LBE模型在ns=0時(shí)，均可簡(jiǎn)化為標(biāo)準(zhǔn)LBE模型，在ns=1時(shí)能夠滿足全反彈LBE模型，能夠滿足質(zhì)量守恒，數(shù)值穩(wěn)定性也較好，因此應(yīng)用比較好。

表1 REV尺度多孔介質(zhì)內(nèi)部分反彈LBE模型（νf為孔隙流體黏度）[11]

1.2 阻力模型

REV尺度多孔介質(zhì)的另一種LBE模型是阻力模型[12-21]，即多孔介質(zhì)的影響通過施加外力得到Stokes或Brinkman方程，這個(gè)思想是基于NITHIARASU等[23]提出的求解不同孔隙度滲流的Navier-Stokes方程，用其他CFD方法（有限差分、有限體積）也可以方便求解，采用LBM方法是眾多CFD方法的一種。該類模型的實(shí)現(xiàn)途徑有兩類：第一類是修正速度[12-17]，即將多孔介質(zhì)的阻礙作用通過修正速度來體現(xiàn)；另一種是流體與多孔介質(zhì)構(gòu)造體-固體的相互作用力直接加入分布函數(shù)的演化方程中[18-21]。

1.2.1 修正速度的阻力模型

修正速度的REV尺度多孔介質(zhì)LBE阻力模型是在碰撞步中修正瞬時(shí)速度來描述多孔介質(zhì)對(duì)流體的影響，所以其演化方程與標(biāo)準(zhǔn)BGK演化方程相同[12]，如式(8)。

1997年SPAID和PHELAN[12]建立了直接求解BRINKMAN方程的LBE模型（SP-LBE模型）。該模型將流體-固體的相互作用力的影響通過修正平衡分布函數(shù)的速度來實(shí)現(xiàn)，這種作用力處理方式對(duì)應(yīng)的宏觀流動(dòng)方程中存在一些誤差，忽略誤差后，SP-LBE對(duì)應(yīng)的宏觀方程在穩(wěn)態(tài)時(shí)就是Brinkman方程。SP-LBE模型平衡態(tài)分布函數(shù)的計(jì)算采用修正后的平衡態(tài)速度，修正的平衡態(tài)速度包含了多孔介質(zhì)的阻力作用，而阻力大小又與當(dāng)?shù)厮俣认嚓P(guān)，如式(9)～式(11)。

式中， β＝v/K，是滲透勢(shì)能參數(shù)，是黏度v與滲透率K的比值。

與SP-LBE模型類同，F(xiàn)REED[13]把多孔介質(zhì)阻力通過碰撞過程在平均速度中體現(xiàn)[Freed-LBE模型，式(13)]，可以用于Mach數(shù)較大的流動(dòng)過程。與多孔介質(zhì)流動(dòng)控制方程相比，F(xiàn)reed-LBE模型對(duì)應(yīng)的宏觀方程存在人工多余項(xiàng)，需要通過具體問題來分析。該模型比較適合用于具有裂隙及導(dǎo)電介質(zhì)的巖石孔隙中流體及其他工程模擬。Freed-LBE模型公式如式(12)～式(14)。

Shan-Chen模型[14-15]最初是應(yīng)用于多相流和多組分流動(dòng)，把體積力F在平衡速度中體現(xiàn)，體積力F主要包括液液流體間相互作用力、液固相間相互作用力和外力。LI等[16-17]基于Shan-Chen模型提出單相單組分LBM方法（Li-LBE模型），其平衡分布函數(shù)可以簡(jiǎn)化成式(16)，平衡態(tài)速度為式(17)。該模型與其他模型相比，不用計(jì)算體積力F，節(jié)省計(jì)算時(shí)間，提高計(jì)算效率。

式中，K為滲透率，ε為孔隙率，G為外部體積力。

1.2.2 修正演化方程的阻力模型

修正演化方程的阻力模型是將多孔介質(zhì)效應(yīng)直接加入到演化方程，其演化方程[19]為式(19)。

式中，F(xiàn)i是F的函數(shù)，根據(jù)其處理方式不同，可以分成不同的阻力模型。

MARTYS[18]對(duì)SP-LBE模型進(jìn)行改進(jìn)（Martys-LBE模型），采用LUO[24]的作用力處理方法描述介質(zhì)阻力，在演化方程中直接增加了一個(gè)作用力項(xiàng)來消除部分離散誤差，但是不能完全消除，如式(20)。其平衡態(tài)分布函數(shù)與SP-LBE模型相同，即式(9)，但是作用力不修正平衡態(tài)分布函數(shù)的速度，而是直接修正演化方程。Martys-LBE模型比SP-LBE模型的數(shù)值穩(wěn)定性好，因?yàn)镾P-LBE模型修正平衡分布函數(shù)的速度，增加了非物理項(xiàng)。Martys-LBE模型其作用力表達(dá)式如式(20)此時(shí)

GUO等[19]在2002年提出了求解包含線性阻力、黏性項(xiàng)和非線性阻力項(xiàng)的通用滲流模型的REV尺度滲流的廣義LBE模型（GLBE模型）。該模型的平衡分布函數(shù)和作用力項(xiàng)中均中包含反映介質(zhì)特性的孔隙率ε，如式(21)和式(22)，而且對(duì)滲流速度大小沒有要求。此模型采用GUO等[25]提出的作用力模型，不存在離散誤差，可得到準(zhǔn)確的宏觀方程。其中，F(xiàn)是包含介質(zhì)阻力和外部體積力的總力，式(23)右端第一和第二項(xiàng)是線性（Darcy）和非線性（Forchheimer）介質(zhì)阻力。Fε是結(jié)構(gòu)參數(shù)，與孔隙率有關(guān)，對(duì)由固體顆粒構(gòu)成的多孔介質(zhì)，如公式(24)。

目前還有一種趨勢(shì)，將LBM與其他數(shù)值方法結(jié)合起來模擬多孔介質(zhì)內(nèi)流動(dòng)和傳遞過程。ZARGHAMI等[20]把GLBE模型和有限體積法相結(jié)合，用有限體積方法離散對(duì)流擴(kuò)散項(xiàng)，如式(25)，并用校正因子來提高數(shù)值穩(wěn)定性。該模型通過Poiseuille流、Couette流和平板驅(qū)動(dòng)流的驗(yàn)證，用于模擬均質(zhì)和非均質(zhì)多孔介質(zhì)管道流。

式中，Ax.y為abcd的面積，ac、ab、bd、bc分別為格子邊界。

ZHANG等[21]提出了Navier-Stokes方程的顯性格式解的LBE模型（Zhang-LBE模型），該模型為通過Chapman-Enskog和泰勒展開，省略高階項(xiàng)之后，得到平均體積Navier-Stokes方程，在演化方程[式(26)]中加入了附加項(xiàng)Pi，如式(28)。該模型只在平衡分布函數(shù)中加入反應(yīng)介質(zhì)特性的孔隙率，如式(27)，但是與GLBE模型中孔隙率加入方式不同，作用力項(xiàng)處理方法采用式(29)，F(xiàn)公式與GLBE模型一致，見式(23)。該模型能有效地解決多孔介質(zhì)流和兩相流。

兩種類型多孔介質(zhì)REV尺度LBE模型，應(yīng)用較為廣泛的是Freed-LBE模型和GLBE模型。如果考慮不可壓流時(shí)，GLBE模型更加精確；對(duì)于穩(wěn)態(tài)不可壓低速多孔介質(zhì)流動(dòng)問題，兩者的區(qū)別在于GLBE模型加入了非線性阻力項(xiàng)，在低速情況下，該項(xiàng)的作用較弱[26]。GLBE模型的作用力是基于GUO等的作用力模型，可以準(zhǔn)確得到宏觀方程，不存在離散誤差；其次，該模型的平衡分布函數(shù)中包含反應(yīng)介質(zhì)特性的孔隙率，而其余REV尺度LBE模型的平衡分布函數(shù)與標(biāo)準(zhǔn)模型類似，不能反映介質(zhì)特性；最后，GLBE模型可以用于高速流動(dòng)。

2　REV尺度LBE模型應(yīng)用現(xiàn)狀

2.1 REV尺度LBE模型應(yīng)用于多孔介質(zhì)內(nèi)流動(dòng)

各REV尺度LBE模型在多孔介質(zhì)內(nèi)流動(dòng)中均有運(yùn)用，但是運(yùn)用最廣泛的是GLBE模型。

WBS-LBE模型、YNS-LBE模型和ZM-LBE模型均被用于模擬不同滲透特性多孔介質(zhì)內(nèi)的流動(dòng)[10-11]，如圖1。模擬結(jié)果顯示W(wǎng)BS-LBE模型在界面處不連續(xù)，YNS-LBE模型和ZM-LBE模型在界面處速度連續(xù)，而且在ns較小時(shí)的連續(xù)性比較好（ns1=ns3=0.2，ns2=0.1），且多孔介質(zhì)的孔隙度變化不會(huì)影響中心處的速度。

圖1 不同滲透特性多孔介質(zhì)示意圖及其速度分布圖[10-11]

REV尺度LBE模型不僅被用于分層多孔介質(zhì)內(nèi)速度的連續(xù)性，還被用于多孔介質(zhì)內(nèi)Klinkenberg效應(yīng)等基礎(chǔ)研究?？紫洞笮≡诩{米到微米范圍內(nèi)，克努森數(shù)（Kn數(shù)，分子平均自由程和流動(dòng)的宏觀特征長(zhǎng)度之比）相對(duì)較大時(shí)，氣體不再滿足連續(xù)介質(zhì)假設(shè)，導(dǎo)致氣體滑移流、瞬變流，甚至真空分子流。由于存在滑移效應(yīng)，測(cè)得表觀滲透率高于固有滲透率，兩者之差隨著Kn數(shù)增大而增大。CHEN等[27]基于四參數(shù)隨機(jī)生成方法（QSGS）并運(yùn)用基本元素建模模型（EBB模型）重新構(gòu)造頁(yè)巖，采用BESKOK和KARNIADAKIS[28]提出的表觀滲透率公式，在Kn數(shù)為常數(shù)時(shí)根據(jù)固有滲透率計(jì)算表觀滲透率，運(yùn)用GLBE模型研究了二維含有Klinkenberg效應(yīng)的多孔介質(zhì)內(nèi)流動(dòng)。高孔隙率（ε=0.8）時(shí)，Klinkenberg效應(yīng)可以忽略，隨著孔隙率減小，Klinkenberg效應(yīng)在低滲透多孔介質(zhì)比較明顯，并且隨著Kn數(shù)增大、Klinkenberg效應(yīng)越明顯，隨著壓力的減小、Klinkenberg效應(yīng)更加明顯（圖2）。隨后CHEN等[29]通過Darcy定律預(yù)測(cè)REV尺度的整個(gè)計(jì)算區(qū)域內(nèi)的有效滲透率，研究多孔介質(zhì)成分、顆粒大小、壓力、滑移及吸附作用對(duì)REV尺度滲透率的影響。在致密儲(chǔ)氣巖層或頁(yè)巖層開采天然氣時(shí)，如果不考慮Klinkenberg效應(yīng)，用表觀滲透率代替固有滲透率，不僅會(huì)導(dǎo)致定量錯(cuò)誤而且會(huì)導(dǎo)致定性錯(cuò)誤。

圖2 含有Klinkenberg效應(yīng)和不含Klinkenberg效應(yīng)的速度分布圖[27]（1psi=6894.76Pa）

除了用于致密巖石中滲流，土壤滲流也是REV尺度多孔介質(zhì)LBE模型的研究熱點(diǎn)。樊火[30]通過繪制滲流場(chǎng)幾何剖面圖或者讀入圖形文件，格點(diǎn)化流場(chǎng)、生成LBM格式的邊界數(shù)據(jù)文件，基于GLBE型開發(fā)了土壤滲流LBM_Seepage軟件平臺(tái)，生成數(shù)據(jù)文件，以直觀可視的等值線或者云圖實(shí)現(xiàn)滲流場(chǎng)速度實(shí)時(shí)可視化，并對(duì)豎向具有不同滲透率的土體進(jìn)行了模擬。結(jié)果表明，最上層土的相對(duì)的滲透能力在滲透初期對(duì)滲流速度的傳播起控制作用，只有經(jīng)過一段時(shí)間后土體才會(huì)表現(xiàn)出等效應(yīng)的滲透能力（圖3）。申林芳等[31]進(jìn)一步運(yùn)用GLBE模型探討了土體在壓力作用下孔隙率、滲透率及滲透壓力等影響因素與滲流速度的相互關(guān)系。

Freed-LBE模型被推廣到非均勻網(wǎng)格中[32]，采用HE等[33]提出的插入格式描述非均勻網(wǎng)格，該模型在流動(dòng)和碰撞過程之間采用插值法，不僅能保證數(shù)值準(zhǔn)確性而且易于實(shí)施。該模型對(duì)多孔介質(zhì)、孔隙和邊界的處理方式一樣，避免了不同尺度采用不同數(shù)學(xué)方法而產(chǎn)生的誤差，與其他多尺度LB模型相比，它消除了內(nèi)部邊界，而且采用阻力模型，因此可以用于大尺度。KANG等[32]在裂隙周圍細(xì)化網(wǎng)格研究裂隙對(duì)流動(dòng)的影響，遠(yuǎn)離裂隙的地方采用粗網(wǎng)格（圖4），模擬了復(fù)雜裂隙和簡(jiǎn)化裂隙的多孔介質(zhì)內(nèi)的流動(dòng)，通過K =預(yù)測(cè)REV尺度多孔介質(zhì)內(nèi)的有效滲透率（表2），并提出可以把該模型用曲線坐標(biāo)和復(fù)合網(wǎng)格的設(shè)想。

圖3 不同滲透率的的邊界示意圖及滲透率K由小到大和K由大到小中心線速度圖[30

REV尺度LBE模型不僅可以模擬含裂隙的均質(zhì)巖石中的流動(dòng)，還能根據(jù)通過高分辨率的掃描圖像確定空腔、黏土、長(zhǎng)石、石英不同成分構(gòu)成的非均質(zhì)砂巖結(jié)構(gòu)，如圖5(a)，基于GLBE模型獲得砂巖樣品內(nèi)的流體流動(dòng)速度，然后根據(jù)達(dá)西定律求得整體的滲透率[34-35]。由圖5(b)可以看出，樣品大于10mm時(shí)整體的滲透率趨于常數(shù)，對(duì)于砂巖樣品，大于10mm時(shí)可以達(dá)到REV尺度的要求。由圖5(c)和圖5(d)可以看出，流場(chǎng)主要受無滲透性的石英和空腔影響，由于各成分的孔隙率不同導(dǎo)致各成分邊界處的流線不閉合。但是CT掃描大孔時(shí)容易忽略狹小的孔喉結(jié)構(gòu)而造成誤差，因而掃描圖像的分辨率必須要高[36]。

圖4 裂隙示意圖、非均勻網(wǎng)格示意圖、計(jì)算結(jié)果[32]

表2 滲透率值[32]（格子單位）

圖5 砂巖樣品

REV尺度多孔介質(zhì)模型除了可以通過滲流速度預(yù)測(cè)多孔介質(zhì)的滲透率外，還可以分析多孔介質(zhì)幾何結(jié)構(gòu)特性及外部條件對(duì)滲流速度的影響。GLBE模型被用于在壓力驅(qū)動(dòng)下的非牛頓流體的電滲流過程[37]，在演化方程中加入了基于Herschel-Bulkley模型的非牛頓流體力項(xiàng)，研究外部電勢(shì)驅(qū)動(dòng)流、離子濃度和電場(chǎng)強(qiáng)度、固體顆粒直徑和孔隙度、冪律指數(shù)和屈服應(yīng)力、德拜長(zhǎng)度對(duì)滲透速度的影響。多孔介質(zhì)內(nèi)的流速隨著冪律指數(shù)和屈服應(yīng)力的減小而增大，壁面處的電勢(shì)主要影響壁面處流速，多孔介質(zhì)的電勢(shì)對(duì)通道中心處流速影響較大。低冪律指數(shù)時(shí)，壁面處的電勢(shì)在對(duì)流速的影響更加敏感，特別是在通道高度與德拜長(zhǎng)度的比值較小時(shí)。

HUANG等[38]對(duì)LBM兩相流模型及其應(yīng)用做了詳細(xì)描述，目前LBM兩相流模型已經(jīng)應(yīng)用于REV尺度多孔介質(zhì)內(nèi)多相流的驅(qū)替過程中。SPAID和PHELAN[39]對(duì)纖維多孔介質(zhì)的流動(dòng)進(jìn)行了計(jì)算，纖維內(nèi)部使用SP-LBE模型，外部使用標(biāo)準(zhǔn)LBE模型，模擬圓形纖維內(nèi)樹脂驅(qū)替空氣多組分流動(dòng)過程。纖維間的流動(dòng)比纖維內(nèi)的流動(dòng)速度快，導(dǎo)致纖維內(nèi)部含有空氣空穴，空穴最初是長(zhǎng)軸垂直于流動(dòng)方向的橢圓形，隨后由于表面張力的影響，從橢圓形變?yōu)閳A形（圖6）。根據(jù)達(dá)西定律預(yù)測(cè)的不飽和滲透率在定性上正確。但是需要更多研究證明纖維內(nèi)的空洞隨時(shí)間消失模擬結(jié)果是否準(zhǔn)確。SCHAAP等[40]應(yīng)用阻力模型中Shan-Chen格式多相流LBE模型，模擬了充滿玻璃珠的多孔介質(zhì)內(nèi)水-空氣及水-石油精體系驅(qū)替實(shí)驗(yàn)。該模型基于顯微層析構(gòu)造的多孔介質(zhì)幾何結(jié)構(gòu)，忽略重力、慣性力和黏性力，只考慮毛細(xì)作用力，設(shè)置合理的表面張力和接觸角，得到壓力與飽和度的關(guān)系。水-空氣體系模擬值與測(cè)量值相差不大，但是水-石油精體系模擬值與測(cè)量值差異很大（圖7），原因是模擬條件與實(shí)驗(yàn)條件不一致，如接觸角。該方法也可以用于微尺度界面現(xiàn)象的模擬。PORTER等[41]也運(yùn)用該LBE模型估算潤(rùn)濕-非潤(rùn)濕相界面的面積，模擬值與實(shí)驗(yàn)值吻合。

圖6 不同時(shí)刻單行圓形纖維樹脂驅(qū)替空氣模擬結(jié)果圖[39]（黑色：空氣；灰色：樹脂）

圖7 壓力與飽和度演化圖[40]

阻力模型中的單松弛GLBE模型推廣到多松弛模型，LIU等[42]運(yùn)用二維D2Q8 MRT-GLBE模型模擬不可壓多孔介質(zhì)流，在低黏性時(shí)的數(shù)值穩(wěn)定性比單松弛模型好。隨后該模型被用于解含有時(shí)間導(dǎo)數(shù)和非線性對(duì)流項(xiàng)平均體積動(dòng)量方程，模擬了在多孔介質(zhì)圓柱靜止時(shí)上下壁面以恒定速度運(yùn)動(dòng)和上下壁面靜止時(shí)管道內(nèi)圓柱多孔介質(zhì)以恒定速度運(yùn)動(dòng)的兩種情況下的流動(dòng)[43]。結(jié)果表明R1和R2的本征相平均速度吻合，二者的相平均速度在流體區(qū)吻合，但是在多孔介質(zhì)區(qū)域內(nèi)有偏差[圖8(b)]，并且隨著孔隙度的減小而增大。多孔介質(zhì)阻力對(duì)本征相平均速度和相平均速度有類似規(guī)律。只有含有本征相平均速度的宏觀方程滿足伽利略不變性。多松弛模型與單松弛模型相比，克服了許多限制條件，比如運(yùn)動(dòng)黏度和體積黏度比不再是定值，但是多松弛模型增加了計(jì)算量。

REV尺度LBM已經(jīng)應(yīng)用到工程滲流問題。何瑩松[44]將REV尺度GLBE模型應(yīng)用于煤礦開采過程中回采工作面瓦斯?jié)B流問題，該研究中認(rèn)為煤體的孔隙度受應(yīng)力重新分布的影響很小，孔隙度基本不發(fā)生變化，滲透率沿流動(dòng)方向呈現(xiàn)分段函數(shù)的不均勻分布，得到了不同時(shí)刻非均質(zhì)煤層中瓦斯壓力分布和速度分布。這些數(shù)據(jù)能夠生動(dòng)再現(xiàn)瓦斯流動(dòng)過程，有助于指導(dǎo)煤礦安全生產(chǎn)。

圖8 多孔介質(zhì)圓柱運(yùn)動(dòng)示意圖及ε=0.7、Re=100通道中心固有相平均速度、相平均速度圖[42]

2.2 REV尺度LBE模型應(yīng)用于多孔介質(zhì)內(nèi)流動(dòng)與傳熱耦合過程

流動(dòng)、傳熱和傳質(zhì)兩個(gè)及兩個(gè)以上過程耦合時(shí)，LBE模型可以分為3種；即多速模型、雙分布函數(shù)模型（DDF）與差分方法相結(jié)合的混合模型。目前文獻(xiàn)大多采用DDF模型來模擬多過程，該模型使用兩個(gè)（及以上）分布函數(shù)演化方程，分別用于速度場(chǎng)和溫度場(chǎng)（濃度場(chǎng)），格子結(jié)構(gòu)簡(jiǎn)單，且具有良好的數(shù)值穩(wěn)定性。

等溫GLBE模型通過DDF模型推廣到多孔介質(zhì)熱流動(dòng)中（DDF-GLBE模型）[45]，速度場(chǎng)用等溫GLBE模型模擬，而溫度場(chǎng)使用一個(gè)新的分布函數(shù)描述。DDF-GLBE模型被用于多孔介質(zhì)方腔模擬孔隙度恒定[45-49]、孔隙度變化[50]、含有熱源[51]的速度場(chǎng)和溫度場(chǎng)，研究孔隙度、雷諾數(shù)（Re數(shù)）、達(dá)西數(shù)（Da數(shù)）等參數(shù)對(duì)速度場(chǎng)和溫度場(chǎng)的影響。DDF-GLBE模型還被用于研究邊界充滿多孔介質(zhì)和中心充滿多孔介質(zhì)兩種情況下管道內(nèi)的流動(dòng)與傳熱過程[52-53]。

REV尺度的GLBE模型可以直接將孔隙尺度模擬得到的滲透率和有效熱導(dǎo)率帶入REV模型，避免了采用經(jīng)驗(yàn)公式可能得到的錯(cuò)誤結(jié)論[54]。模擬結(jié)果表明，REV尺度得到的熱壁上的平均Nusselt數(shù)與孔隙尺度下的趨勢(shì)一致，證明REV尺度在定性上是正確的。在瑞利數(shù)（Ra數(shù)）低于107時(shí)，GLBE模型能定量預(yù)測(cè)孔隙尺度下自然對(duì)流的宏觀規(guī)律，Ra數(shù)大的時(shí)候只能定性預(yù)測(cè)孔隙尺度的宏觀規(guī)律（圖9）。REV尺度GLBE模型對(duì)滲透率如何影響多孔介質(zhì)自然對(duì)流宏觀現(xiàn)象刻畫很好，與孔隙尺度模擬結(jié)果一致，但是往往夸大了流固熱導(dǎo)率帶來的影響。

圖9 不同Ra數(shù)下熱壁上平均Nusselt數(shù)隨方柱個(gè)數(shù)N的變化結(jié)果對(duì)比[54]

REV尺度的LBE模型為太陽(yáng)墻工程應(yīng)用研究提供參考，通過含熱源項(xiàng)的LBE模型對(duì)太陽(yáng)墻系統(tǒng)中復(fù)雜流場(chǎng)及熱傳遞進(jìn)行模擬[55-56]，可以預(yù)測(cè)太陽(yáng)墻在不同工作條件下的狀況，并優(yōu)化太陽(yáng)墻系統(tǒng)的設(shè)計(jì)。

針對(duì)工程中常見的多孔介質(zhì)圓管內(nèi)傳熱問題，榮伏梅等[57-58]構(gòu)造了REV尺度軸對(duì)稱熱DDF-GLBE模型，包含弱可壓軸對(duì)稱LBE模型和不可壓軸對(duì)稱LBE模型，把刻畫多孔介質(zhì)的參數(shù)加入到分布函數(shù)中并構(gòu)造合適的分布函數(shù)，使之能恢復(fù)正確的宏觀方程，具有形式簡(jiǎn)單、外力項(xiàng)不含復(fù)雜梯度形式的優(yōu)點(diǎn)，并應(yīng)用于非線性滲流數(shù)值模擬。分別研究了覆蓋多孔介質(zhì)的方柱繞流、充滿多孔介質(zhì)內(nèi)的變黏性Rayleigh-Benard對(duì)流、以及圓管內(nèi)插入多孔介質(zhì)強(qiáng)化對(duì)流傳熱非線性滲流問題。由于這類問題的計(jì)算量都很大，單個(gè)CPU的計(jì)算效率較低，因此采用GPU（圖形處理器）加速計(jì)算，可支持復(fù)雜的計(jì)算任務(wù)[59]。

單松弛DDF-GLBE模型推廣到多松弛模型中，建立REV尺度二維多松弛（MRT）DDF-GLBE方程，用D2Q9 MRT-GLBE模型模擬流場(chǎng)，同時(shí)用D2Q5 MRT-LBE模型模擬溫度場(chǎng)，模擬結(jié)果和解析解與文獻(xiàn)計(jì)算結(jié)果吻合[60]，而且在低黏度時(shí)，MRT-LBE模型的數(shù)值穩(wěn)定性比LBGK模型好，但是計(jì)算量大。該直角坐標(biāo)系下的MRT-DDF-GLBE模型延伸到軸對(duì)稱模型[61]。

2.3 REV尺度LBE模型應(yīng)用于多孔介質(zhì)內(nèi)流動(dòng)、傳質(zhì)和化學(xué)反應(yīng)耦合過程

與流動(dòng)和傳熱過程類同，仍使用DDF模型描述流動(dòng)和傳質(zhì)兩個(gè)過程，即運(yùn)用兩個(gè)分布函數(shù)分別模擬流場(chǎng)和濃度場(chǎng)，并在濃度演化方程中加入反應(yīng)源項(xiàng)，用于模擬化學(xué)反應(yīng)。TIAN等[62]基于DDF-GLBE模型模擬了CO2地下填埋地球化學(xué)反應(yīng)過程，研究CO2注入過程中溶質(zhì)離子反應(yīng)及擴(kuò)散、巖石骨架CaCO3溶解相互作用，引起孔隙率和滲透率變化，反作用于流場(chǎng)。如圖10，開始在一段時(shí)間內(nèi)出口處的H+、Ca2+和HCO3?離子濃度為0，隨著反應(yīng)的進(jìn)行，3種離子濃度逐漸增加。由于反應(yīng)消耗和對(duì)流擴(kuò)散的相互作用，在出口處反應(yīng)物H+離子濃度在1×105時(shí)出現(xiàn)變化。Ca2+離子濃度到達(dá)最大后保持平衡，是因?yàn)镃a2+離子是可逆反應(yīng)。

該模型還用于研究有機(jī)廢水溶液中化學(xué)反應(yīng)、燃料電池和質(zhì)子交換膜燃料電池中[56，63-66]。LIAO 等[63]和楊艷霞[64]研究了有機(jī)廢水中含擾流光合細(xì)菌包埋顆粒的流動(dòng)、傳質(zhì)及其內(nèi)部的光生化反應(yīng)過程，分析光照強(qiáng)度、進(jìn)口流速、包埋顆粒的滲透率及孔隙率對(duì)流場(chǎng)、濃度場(chǎng)及產(chǎn)氫性能的影響。XU 等[65-66]把該模型應(yīng)用于固體燃料電池多組分化學(xué)反應(yīng)，計(jì)算濃差極化電勢(shì)，并研究了電池電極厚度、電解質(zhì)濃度、孔隙率、燃料組成等結(jié)構(gòu)參數(shù)和工作溫度操作參數(shù)對(duì)濃差極化、超電勢(shì)、電化學(xué)效率、電池電壓和功率密度的影響，LBM模擬結(jié)果與其他模擬結(jié)果相比更接近實(shí)驗(yàn)值。該模型也被用于模擬質(zhì)子交換膜燃料電池，可以得到電池內(nèi)部各組分濃度的分布狀況[56]。

除了通過DDF模型耦合流動(dòng)和傳質(zhì)外，REV尺度LBE模型耦合偏微分方程和常微分方程的運(yùn)輸方程，模擬了多孔介質(zhì)內(nèi)的擴(kuò)散反應(yīng)過程，得到飽和度和濃度之間的關(guān)系[67]。

REV尺度多分布LBE模型還可以通過DDF模型實(shí)現(xiàn)流動(dòng)、傳熱和傳質(zhì)三者耦合，應(yīng)用到生物過濾池的非等溫過程[68]。比較3種不同結(jié)構(gòu)的生物過濾器，模擬結(jié)果發(fā)現(xiàn)Ra數(shù)對(duì)生物過濾池的凈化率影響很大，調(diào)整Ra數(shù)可以優(yōu)化生物過濾器的凈化效率。停留時(shí)間的標(biāo)準(zhǔn)偏差較小時(shí)生物過濾器內(nèi)的流動(dòng)更加均勻，當(dāng)Ra數(shù)小于臨界值時(shí)，生物過濾器不均質(zhì)性隨著Ra數(shù)的增大而減小，否則，生物過濾器的不均質(zhì)隨著Ra數(shù)的增大而增大（圖11）。

2.4 REV尺度LBE模型應(yīng)用于多孔介質(zhì)相變過程

圖10 出口離子濃度圖[62]

圖11 3種生物過濾器結(jié)構(gòu)圖和停留時(shí)間隨Ra數(shù)變化的標(biāo)準(zhǔn)平均誤差[68]

REV尺度LBE模型通過加入相變熱對(duì)溫度的影響，可以處理相變過程。流場(chǎng)基于GLBE模型，而溫度演化方程中加相變焓模型，可以用于凝固、冷凝、和融化等過程。該模型被用于模擬多孔介質(zhì)內(nèi)自然對(duì)流對(duì)融化的影響過程[69]，LBM的模擬結(jié)果和前人的實(shí)驗(yàn)及模擬結(jié)果吻合（圖12b），Ra數(shù)、Da數(shù)和孔隙率增大時(shí)，自然對(duì)流對(duì)融化過程的影響作用也增強(qiáng)。該模型還用于模擬固化及融化過程中Ra數(shù)對(duì)溫度和流動(dòng)的影響及孔隙度對(duì)對(duì)流傳熱和相界面位置的影響[55，70]。類似模型用于探究了密封的同心圓環(huán)內(nèi)充滿Al2O3的多孔介質(zhì)內(nèi)冰的融化過程，該模型把相變熱作為源項(xiàng)加到溫度演化方程中[71]。

圖12 多孔介質(zhì)內(nèi)從左向右融化示意圖和不同時(shí)刻界面圖[69]

REV尺度LBE模型用于模擬了各向同性的均質(zhì)充滿水的可變形多孔介質(zhì)內(nèi)的蒸發(fā)干燥過程[72]，該模型分別用不同的演化方程處理液相、固相變形速度、壓力及溫度場(chǎng)，分析了泊松系數(shù)、楊氏模量和滲透率對(duì)宏觀流場(chǎng)的影響。

REV尺度單松弛的相變模型推廣到多松弛相變模型，模擬多孔介質(zhì)內(nèi)瞬態(tài)固液相變過程[73]，該模型用Brinkman-Forchheimer的達(dá)西擴(kuò)展模型模擬動(dòng)量傳遞，把相變焓作為源項(xiàng)加到MRT-LBM的溫度演化方程中來確定固液的相界面。與BGK-LBE模型相比，在低黏度時(shí)數(shù)值穩(wěn)定性好，并且邊界條件處理更精確。

3　結(jié)論和展望

化工過程、油藏氣開發(fā)、環(huán)境保護(hù)及新能源開發(fā)都廣泛涉及多孔介質(zhì)內(nèi)的流動(dòng)、傳遞和反應(yīng)過程。REV尺度LBE模型無需多孔介質(zhì)內(nèi)的具體結(jié)構(gòu)、且計(jì)算域大，比較適合解決包含多孔介質(zhì)結(jié)構(gòu)的地下水滲流、污染物遷移、固定床反應(yīng)器等工程應(yīng)用問題。在眾多的部分反彈模型和阻力模型中，由于GLBE模型的平衡分布函數(shù)和作用力項(xiàng)中均包含反映介質(zhì)特性的孔隙度，且不含離散誤差，因而應(yīng)用廣泛。直角坐標(biāo)系下的REV尺度多孔介質(zhì)的單松弛LBE模型已經(jīng)推廣到了多松弛模型和軸對(duì)稱模型，REV尺度多孔介質(zhì)的多松弛模型可以用于低粘度流體，有更多可調(diào)參數(shù)，且松弛因子與宏觀方程無關(guān)，計(jì)算更穩(wěn)定，但計(jì)算量要比單松弛模型增加15%左右。

多孔介質(zhì)REV尺度的LBE模型可以用于模擬均質(zhì)、分層均質(zhì)和含有裂隙的多孔介質(zhì)內(nèi)的流動(dòng)，也可以與掃描電鏡等構(gòu)造多孔介質(zhì)方法相結(jié)合用于模擬非均質(zhì)的多孔介質(zhì)內(nèi)的流動(dòng)及電滲流等，或者與QSGS方法和EBB方法相結(jié)合研究Klinkenberg效應(yīng)。應(yīng)用多分布LBE模型可探索多孔介質(zhì)內(nèi)的流動(dòng)、傳熱、傳質(zhì)及化學(xué)反應(yīng)等多過程問題，還可以用于融化、凝固和蒸發(fā)等相變過程。REV尺度多孔介質(zhì)的LBE模型已經(jīng)成功用于非牛頓流體的電滲流過程、固定床反應(yīng)器多孔介質(zhì)內(nèi)多相流動(dòng)、煤層氣瓦斯?jié)B漏、土壤滲流、太陽(yáng)墻系統(tǒng)設(shè)計(jì)、石油巖層多相驅(qū)替和燃料電池等工程問題中。

目前REV尺度多孔介質(zhì)LBE模型模擬主要集中在二維尺度，為了更加接近實(shí)際物理過程，今后應(yīng)該把二維REV多孔介質(zhì)LBE模型推廣到三維模型中，把直角坐標(biāo)系的LBE模型延展到曲線坐標(biāo)，更加符合實(shí)際多孔介質(zhì)。對(duì)于部分反彈模型可以把壓力梯度引入到固體散射密度中，使多孔介質(zhì)的滲透特性隨著壓力梯度的變化而變化，用于模擬非達(dá)西流。

但是多孔介質(zhì)內(nèi)發(fā)生的三傳一反過程，在更小的介觀尺度上，多孔介質(zhì)幾何結(jié)構(gòu)是非均質(zhì)的[74-75]，在更大工業(yè)尺度上，多孔介質(zhì)內(nèi)物性參數(shù)及過程呈現(xiàn)出各向異性，為了更科學(xué)準(zhǔn)確地描述多孔介質(zhì)內(nèi)發(fā)生的各個(gè)過程，不能忽略孔隙尺度參數(shù)對(duì)滲透率、傳熱傳質(zhì)系數(shù)的影響，更要考慮更大工程尺度上各個(gè)過程各向異性的觸發(fā)機(jī)理和演化過程。根據(jù)作者課題組所完成的工作，期望在孔隙尺度上獲得包含個(gè)孔隙幾何結(jié)構(gòu)參數(shù)的滲透率、有效傳熱、傳質(zhì)系數(shù)，將其帶入到REV尺度多孔介質(zhì)的LBE數(shù)學(xué)模型中，可避免經(jīng)驗(yàn)公式帶來的誤差. 同時(shí)從LBM演化方程本身出發(fā)，建立一個(gè)能夠具有明確合理物理意義的各向異性演化方程，能夠獲得更大尺度上包含各向異性特性的溫度和濃度的分布特性。因而建立起基于LBM多尺度數(shù)學(xué)模型和數(shù)值方法及其求解平臺(tái)。

符 號(hào) 說 明

Ax.y—— abcd的面積

Da—— 達(dá)西數(shù)

eq —— 平衡態(tài)

F—— 體積力

Fε—— 結(jié)構(gòu)參數(shù)

fi—— 分布函數(shù)

fi—— 反分布函數(shù)

fic—— 碰撞后的分布函數(shù)

G—— 外部體積力

i—— 速度離散方向

K—— 滲透率

ns—— 固體格點(diǎn)的平均散射密度

Pi—— 附加項(xiàng)

Ra—— 瑞利數(shù)

Re—— 雷諾數(shù)

t—— 離散時(shí)間

t*—— 碰撞前的時(shí)間

u'—— 碰撞后速度

u(eq)—— 平衡態(tài)速度

νf—— 孔隙流體黏度

Δi(x,t)—— 流體粒子與固體骨架之間的碰撞

β —— 滲透勢(shì)能參數(shù)

ν —— 黏度

ρ—— 密度

ε—— 孔隙率

τ—— 松弛時(shí)間

ω—— 權(quán)系數(shù)

參 考 文 獻(xiàn)

[1] 何雅玲，王勇，李慶. 格子Boltzmann方法的理論及應(yīng)用 [M]. 北京：科學(xué)出版社，2009.

[2] SUCCI S，BENZI R，MASSAIOLI F. A review of the lattice Boltzmann method [J]. International Journal of Modern Physics C，1993，4（2）：409-415.

[3] CHEN S Y，DOOLEN G D. Lattice Boltzmann method for fluid flows [J]. Annual Review of Fluid Mechanics，1998，30（1）：329-364.

[4] 李楠，史俊瑞，肖紅俠，等. 格子Boltzmann方法研究及其應(yīng)用 [J].沈陽(yáng)工程學(xué)院學(xué)報(bào)(自然科學(xué)版)，2015（2）：116-118，122.

[5] SUKOP M C，THORNE D T. Lattice Boltzmann modeling：an introduction for geoscientists and engineers [M]. Springer，2006.

[6] GAO Y，SHARMA M. A LGA model for fluid flow in heterogeneous porous media [J]. Transport in Porous Media，1994，17（1）：1-17.

[7] DARDIS O，MCCLOSKEY J. Lattice Boltzmann scheme with real numbered solid density for the simulation of flow in porous media [J]. Physical Review E，1998，57（4）：4834-4837.

[8] THORNE D T，SUKOP M C. Lattice Boltzmann model for the Elder problem [J]. Developments in Water Science，2004，55：1549-1557.

[9] WALSH S D，BURWINKLE H，SAAR M O. A new partial-bounceback lattice-Boltzmann method for fluid flow through heterogeneous media [J]. Computers & Geosciences，2009，35（6）：1186-1193.

[10] YEHYA A，NAJI H，SUKOP M. Simulating flows in multi-layered and spatially-variable permeability media via a new Gray Lattice Boltzmann model [J]. Computers and Geotechnics，2015，70：150-158.

[11] ZHU J J，MA J S. An improved gray lattice Boltzmann model for simulating fluid flow in multi-scale porous media [J]. Advances in Water Resources，2013，56：61-76.

[12] SPAID M A，PHELAN JR F R. Lattice Boltzmann methods for modeling microscale flow in fibrous porous media [J]. Physics of Fluids，1997，9（9）：2468-2474.

[13] FREED D M. Lattice-Boltzmann method for macroscopic porous media modeling [J]. International Journal of Modern Physics C，1998，9（8）：1491-1503.

[14] KANG Q，ZHANG D，CHEN S. Displacement of a two-dimensional immiscible droplet in a channel [J]. Physics of Fluids，2002，14（9）：3203-3214.

[15] SHAN X W，CHEN H D. Lattice Boltzmann model for simulating flows with multiple phases and components [J]. Physical Review E，1993，47（3）：1815.

[16] LI J，BROWN D. Upscaled Lattice Boltzmann method for simulations of flows in heterogeneous porous media [J]. arXiv preprint arXiv：13120195，2013，

[17] ILIEV O，EFENDIEV Y，LI J，et al. Multiscale Lattice Boltzmann Method for flow simulations in highly heterogenous porous media[C]//proceedings of the SPE Reservoir Characterization and Simulation Conference and Exhibition，Society of Petroleum Engineers，2013.

[18] MARTYS N S. Improved approximation of the Brinkman equation using a lattice Boltzmann method [J]. Physics of Fluids，2001，13（6）：1807-1810.

[19] GUO Z L，ZHAO T S. Lattice Boltzmann model for incompressible flows through porous media [J]. Physical Review E，2002，66（3）：036304.

[20] ZARGHAMI A，BISCARINI C，SUCCI S，et al. Hydrodynamics in porous media：a finite volume lattice Boltzmann study [J]. Journal of Scientific Computing，2014，59（1）：80-103.

[21] ZHANG J F，WANG L M，OUYANG J. Lattice Boltzmann model for the volume-averaged Navier-Stokes equations [J]. Europhysics Letters，2014，107（3）：20001.

[22] CHEN Y L，ZHU K Q. A study of the upper limit of solid scattersdensity for gray Lattice Boltzmann Method [J]. Acta Mechanica Sinica，2008，24（5）：515-522.

[23] NITHIARASU P，SEETHARAMU K，SUNDARARAJAN T. Natural convective heat transfer in a fluid saturated variable porosity medium [J]. International Journal of Heat and Mass Transfer，1997，40（16）：3955-3967.

[24] LUO L S. Unified theory of lattice Boltzmann models for nonideal gases [J]. Physical Review Letters，1998，81（8）：1618.

[25] GUO Z L，ZHENG C G，SHI B C. Discrete lattice effects on the forcing term in the lattice Boltzmann method [J]. Physical Review E，Statistical，Nonlinear，and Soft Matter Physics，2002，65(4 Pt 2B)：046308.

[26] 雷鳴. 關(guān)于兩種REV尺度多孔介質(zhì)LBE模型的討論 [J]. 科學(xué)技術(shù)與工程，2011（20）：4814-4820.

[27] CHEN L，F(xiàn)ANG W Z，KANG Q J，et al. Generalized lattice Boltzmann model for flow through tight porous media with Klinkenberg's effect [J]. Physical Review E，2015，91（3）：033004.

[28] BESKOK A，KARNIADAKIS G E. Report：a model for flows in channels，pipes，and ducts at micro and nano scales [J]. Microscale Thermophysical Engineering，1999，3（1）：43-77.

[29] CHEN L，KANG Q，DAI Z，et al. Permeability prediction of shale matrix reconstructed using the elementary building block model [J]. Fuel，2015，160：346-356.

[30] 樊火. 土體飽和滲流的Lattice Boltzmann Method數(shù)值模擬研究[D]. 長(zhǎng)沙：長(zhǎng)沙理工大學(xué)，2009.

[31] 申林方，王志良，李邵軍. 基于格子博爾茲曼方法表征體元尺度土體細(xì)觀滲流場(chǎng)的數(shù)值模擬 [J]. 巖土力學(xué)，2015（s2）：689-694.

[32] KANG Q J，ZHANG D X，CHEN S Y. Unified lattice Boltzmann method for flow in multiscale porous media [J]. Physical Review E，2002，66(5)：056307.

[33] HE X Y，LUO L S，DEMBO M. Some progress in lattice Boltzmann method. Part I. Nonuniform mesh grids [J]. Journal of Computational Physics，1996，129（2）：357-363.

[34] GAO J F，XING H L，LI Q，et al. Scale effect of 3d heterogeneous porous media on geo-fluid characteristics[C]//Insight from Massively Parallel Lattice Boltzmann Computing. proceedings of the SPE Unconventional Resources Conference and Exhibition-Asia Pacific，Society of Petroleum Engineers ，2013.

[35] GAO J，XING H，TIAN Z，et al. Lattice Boltzmann modeling and evaluation of fluid flow in heterogeneous porous media involving multiple matrix constituents [J]. Computers & Geosciences，2014，62：198-207.

[36] GAO J，XING H，RUDOLPH V，et al. Parallel lattice boltzmann computing and applications in core sample feature evaluation [J]. Transport in Porous Media，2015，107（1）：65-77.

[37] TANG G H，YE P X，TAO W Q. Pressure-driven and electroosmotic non-Newtonian flows through microporous media via lattice Boltzmann method [J]. Journal of Non-Newtonian Fluid Mechanics，2010，165（21）：1536-1542.

[38] HUANG H，SUKOP M，LU X. Multiphase lattice Boltzmann methods：theory and application [M]. New York：John Wiley & Sons，2015.

[39] SPAID M A A，PHELAN F R. Modeling void formation dynamics in fibrous porous media with the lattice Boltzmann method [J]. Compos Part A-Appl S，1998，29（7）：749-755.

[40] SCHAAP M G，PORTER M L，CHRISTENSEN B S，et al. Comparison of pressure‐saturation characteristics derived from computed tomography and lattice Boltzmann simulations [J]. Water Resources Research，2007，43（12）.DOI：10.1029/2006WR005730.

[41] PORTER M L，SCHAAP M G，WILDENSCHILD D. Lattice-Boltzmann simulations of the capillary pressure–saturation–interfacial area relationship for porous media [J]. Advances in Water Resources，2009，32(11)：1632-1640.

[42] LIU Q，HE Y L. Multiple-relaxation-time lattice Boltzmann modeling of incompressible flows in porous media [J]. Physica A：Statistical Mechanics and its Applications，2015，429：215-230.

[43] WANG L，WANG L P，GUO Z，et al. Volume-averaged macroscopic equation for fluid flow in moving porous media [J]. International Journal of Heat and Mass Transfer，2015，82：357-368.

[44] 何瑩松. 基于格子Boltzmann方法的多孔介質(zhì)流體滲流模擬 [J].科技通報(bào)，2013（4）：118-120.

[45] GUO Z L，ZHAO T S. A lattice Boltzmann model for convection heat transfer in porous media [J]. Numerical Heat Transfer，Part B，2005，47(2)：157-177.

[46] XU Y S，WANG R M，LI G N，et al. Numerical simulation of driven convective heat transfer based lattice boltzmann method in a porous cavity [J]. Journal of Chemistry，2015，501：710341.

[47] SETA T，TAKEGOSHI E，OKUI K. Lattice Boltzmann simulation of natural convection in porous media [J]. Mathematics and Computers in Simulation，2006，72（2）：195-200.

[48] CHAI Z H，GUO Z L，SHI B C. Lattice Boltzmann simulation of mixed convection in a driven cavity packed with porous medium [C]//computational science–ICCS 2007，Springer，2007：802-809.

[49] 萬斌，羅成龍. 格子-Boltzmann方法模擬多孔介質(zhì)內(nèi)自然對(duì)流蓄熱過程 [J]. 南昌大學(xué)學(xué)報(bào)(工科版)，2012（2）：165-167+171.

[50] YAN W W，LIU Y，GUO Z L，et al. Lattice Boltzmann simulation on natural convection heat transfer in a two-dimensional cavity filled with heterogeneously porous medium [J]. International Journal of Modern Physics C，2006，17（6）：771-783.

[51] 許超敏. 用格子Boltzmann方法模擬熱流動(dòng) [D]. 長(zhǎng)春：吉林大學(xué)，2010.

[52] ZARGHAMI A，F(xiàn)RANCESCO S D，BISCARINI C. Porous substrate effects on thermal flows through a REV-scale finite volume lattice Boltzmann model [J]. International Journal of Modern Physics C，2014，25（2）：1350086.

[53] SHOKOUHMAND H，JAM F，SALIMPOUR M. Simulation of laminar flow and convective heat transfer in conduits filled with porous media using Lattice Boltzmann Method [J]. International Communications in Heat and Mass Transfer，2009，36（4）：378-384.

[54] 魯建華. 基于格子Boltzmann方法的多孔介質(zhì)內(nèi)流動(dòng)與傳熱的微觀模擬 [D]. 武漢：華中科技大學(xué)，2009.

[55] 邵九姑. 用格子Boltzmann方法模擬研究多孔介質(zhì)中的傳質(zhì)傳熱問題 [D]. 杭州：浙江師范大學(xué)，2012.

[56] 陳柔. 格子Boltzmann方法的若干應(yīng)用 [D]. 杭州：浙江師范大學(xué)，2013.

[57] RONG F M，GUO Z L，CHAI Z H，et al. A lattice Boltzmann model for axisymmetric thermal flows through porous media [J]. International Journal of Heat and Mass Transfer，2010，53（23）：5519-5527.

[58] 榮伏梅. 基于格子Boltzmann方法的多孔介質(zhì)REV尺度的流動(dòng)與強(qiáng)化傳熱研究 [D]. 武漢：華中科技大學(xué)，2011.

[59] RONG F M，ZHANG W H，SHI B C，et al. Numerical study of heat transfer enhancement in a pipe filled with porous media by axisymmetric TLB model based on GPU [J]. International Journal of Heat and Mass Transfer，2014，70：1040-1049.

[60] LIU Q，HE Y L，LI Q，et al. A multiple-relaxation-time lattice Boltzmann model for convection heat transfer in porous media [J]. International Journal of Heat and Mass transfer，2014，73：761-775.

[61] LIU Q，HE Y L，LI Q. Multiple-relaxation-time lattice Boltzmann model for simulating axisymmetric thermal flows in porous media [J]. arXiv preprint arXiv：14095929，2014.

[62] TIAN Z W，XING H L，TAN Y L，et al. A coupled lattice Boltzmann model for simulating reactive transport in CO2injection [J]. Physica A：Statistical Mechanics and its Applications，2014，403：155-164.

[63] LIAO Q，YANG Y X，ZHU X，et al. Lattice Boltzmann simulation of substrate solution through a porous granule immobilized PSB-cell for biohydrogen production [J]. International Journal of Hydrogen Energy，2013，38（35）：15700-15709.

[64] 楊艷霞. 含光生化反應(yīng)的多孔介質(zhì)內(nèi)流動(dòng)及傳輸特性的格子Boltzmann模擬 [D]. 重慶：重慶大學(xué)，2013.

[65] XU H，DANG Z，BAI B F. Numerical simulation of multispecies mass transfer in a SOFC electrodes layer using lattice Boltzmann method [J]. Journal of Fuel Cell Science and Technology，2012，9 （6）：061004.

[66] XU H，DANG Z，BAI B F. Electrochemical performance study of solid oxide fuel cell using lattice Boltzmann method [J]. Energy，2014，67：575-583.

[67] XU Y S，ZHONG Y J，HUANG G X. Lattice Boltzmann method for diffusion-reaction-transport processes in heterogeneous porous media [J]. Chinese Physics Letters，2004，21（7）：1298.

[68] YAN W W，SU Z D，ZHANG H J. Effect of non‐isothermal condition on heterogeneous flow through biofilter media by lattice Boltzmann simulation [J]. Journal of Chemical Technology and Biotechnology，2013，88（3）：456-461.

[69] GAO D Y，CHEN Z Q. Lattice Boltzmann simulation of natural convection dominated melting in a rectangular cavity filled with porous media [J]. International Journal of Thermal Sciences，2011，50（4）：493-501.

[70] 李霞. 用Lattice Boltzmann方法處理傳熱傳質(zhì)中的固液相變問題[D]. 杭州：浙江師范大學(xué)，2011.

[71] JOURABIAN M，F(xiàn)ARHADI M，DARZI A A R. Lattice Boltzmann investigation for enhancing the thermal conductivity of ice using Al2O3porous matrix [J]. International Journal of Computational Fluid Dynamics，2012，26（9/10）：451-462.

[72] EL ABRACH H，DHAHRI H，MHIMID A. Numerical Simulation of Drying of a Saturated Deformable Porous Media by the Lattice Boltzmann Method [J]. Transport in Porous Media，2013，99（3）：427-452.

[73] LIU Q，HE Y L. Double multiple-relaxation-time lattice Boltzmann model for solid-liquid phase change with natural convection in porous media [J]. arXiv preprint arXiv：150308406，2015.

[74] 蔡睿賢，茍晨華，張娜. 多孔介質(zhì)中考慮各向異性自然對(duì)流Brinkman模型的解析解(Ⅱ) [J]. 中國(guó)科學(xué)G輯：物理學(xué)、力學(xué)、天文學(xué)，2005（2）：113-120.

[75] 蔡睿賢，張娜. 多孔介質(zhì)中自然對(duì)流Brinkman模型的解析解 [J].中國(guó)科學(xué)（A輯），2002（6）：566-573.

Models and application of lattice Boltzmann method at REV-scale in porous media

ZHANG Xiaodan1，2，YONG Yumei2，LI Wenjun3，ZHAO Yuansheng4，LI Yuanyuan2，YANG Qiaowen1，YANG Chao2
（1School of Chemical & Environmental Engineering，China University of Mining & Technology (Beijing)，Beijing 100083，China；2Key Laboratory of Green Process and Engineering，Institute of Process Engineering，Chinese Academy of Sciences，Beijing 100190，China；3School of Environmental Engineering，North China Institute of Science and Technology，Langfang 065201，Hebei，China；4Laboratory of Residue Hydrotreating，Research Institute of Petroleum Processing，PetroChina，Beijing 102200，China）

Abstract：This paper discusses the lattice Boltzmann model at representative elementary volume (REV) scale for porous media. According to different treatments of porous media，the lattice Boltzmann model at REV-scale for porous media can be classified into two categories，the partially bouncing-back model and the resistance model. The advantages and disadvantages of various models are analyzed. The Generalized lattice Boltzmann equation (GLBE model) in the resistance model is most widely used. Firstly，the force item of the GLBM model is based on the method proposed by Guo et al，which can beaccurately recovered to the macroscopic equation without discretization error. Secondly，the equilibrium distribution function and the force items involve porosity，which reflects the characteristics of porous media. This review introduces the application of the lattice Boltzmann model at REV-scale in porous media with flow，heat transfer，mass transfer，chemical reaction and phase transition. We should take the pore scale factors into consideration when studying momentum transport，heat transport，mass transport and reaction engineering in porous media at REV-scale. The anisotropic nature of process parameters should be considered at the larger engineering scale. It gives some predictions and perspectives of applications for the lattice Boltzmann model at REV-scale.

Key words：porous media；representative elementary volume (REV)；lattice Boltzmann method (LBM)；flow；heat transfer；mass transfer

中圖分類號(hào)：TQ021.9

文獻(xiàn)標(biāo)志碼：A

文章編號(hào)：1000–6613（2016）06–1698–15

DOI：10.16085/j.issn.1000-6613.2016.06.010

收稿日期：2016-01-19；修改稿日期：2016-01-31。

基金項(xiàng)目：國(guó)家自然科學(xué)基金（21276256,21490584）、國(guó)家重點(diǎn)基礎(chǔ)研究發(fā)展計(jì)劃（2012CB224806）及中央高校基本科研業(yè)務(wù)費(fèi)（3142013097）項(xiàng)目。