格子Boltzmann数值计算的性能优化

格子波兹曼方法
格子波兹曼方法（Lattice Boltzmann Method, LBM）是一种广泛应用于计算流
体力学领域的数值方法。

它基于分子动力学模型，通过离散化空间网格和时间步长来模拟复杂的流体流动问题。

格子波兹曼方法通过将流体宏观物理量离散化到网格上的节点，使用分布函数
描述流体粒子的运动。

流体粒子在相邻节点之间以一种特定的方式进行碰撞和传播，模拟流体的宏观行为。

格子波兹曼方法相对于传统的Navier-Stokes方程求解方法具有多个优势。

首先，它因其并行化的能力而广泛应用于高性能计算中。

其次，LBM的离散化框架使得
它在处理具有复杂边界条件和多相流问题时更加灵活。

此外，LBM对于非连续和
非均匀流体介质的模拟效果也相对较好。

格子波兹曼方法在各个领域都有广泛的应用。

在流体力学领域，LBM被用于
模拟自由表面流动、湍流现象和多孔介质中的流动行为。

在微观领域，LBM也被
用于模拟微观流体力学现象，例如微管流动和纳米颗粒悬浮体的输运行为。

除了流体力学领域，格子波兹曼方法还被应用于其他科学领域。

例如，它被用
于模拟热传导、传质过程、相变以及复杂物质的输运现象。

此外，LBM还被用于
模拟生物流体力学、地下水流动、大气动力学和地震波传播等问题。

综上所述，格子波兹曼方法是一个高效且灵活的数值方法，用于模拟复杂的流
体流动问题。

它在计算流体力学领域以及其他科学领域都有广泛的应用前景。

这种方法的进一步发展和应用将有助于我们更好地理解和预测流体行为，并解决相关领域的实际问题。

颗粒群绕流传质特性的格子Boltzmann模拟王利民;吴承优【摘要】为考察气固两相流中颗粒团聚对气固传质过程的影响,采用耦合传质的格子Boltzmann模型对在空间均匀分布的颗粒结构和非均匀分布的颗粒团聚物结构条件下的绕流过程进行了直接数值模拟,得到了气体绕颗粒流动的速度分布和伴随流动的浓度分布,发现颗粒团聚物对绕流过程的速度分布和浓度分布都具有明显的影响.通过对颗粒绕流的气固传质过程进行定量分析发现,在两种结构条件下计算得到的传质舍伍德数均随着雷诺数的增大而呈指数函数形式增大,并且在均匀结构条件下的传质舍伍德数一般为非均匀结构条件下的3～5倍,颗粒团聚物的存在将严重影响颗粒绕流过程的气固传质效率.通过格子Boltzmann方法建立的气固相间传质模型,可以为研究气固两相流介尺度结构的传递特性提供理论依据.【期刊名称】《化学反应工程与工艺》【年(卷),期】2016(032)004【总页数】8页(P289-296)【关键词】颗粒群绕流;颗粒团聚物;格子玻尔兹曼方法;传质;数值模拟【作者】王利民;吴承优【作者单位】中国科学院过程工程研究所多相复杂系统国家重点实验室,北京100190;中国科学院过程工程研究所多相复杂系统国家重点实验室,北京100190;江西理工大学资源与环境工程学院,江西赣州341000【正文语种】中文【中图分类】TQ021.4气固流动中，团聚物对气固传递特性的影响很大［1］。

从20世纪80年代，李静海等［2］开始研究气固反应器层次的介尺度问题，从颗粒聚团现象入手，认为聚团的形成来源于气体和颗粒的各自运动趋势在竞争中的协调，从而建立了能量最小多尺度（Energy-Minimization Multi-Scale，EMMS）模型，通过EMMS模型计算发现稀相和密相的曳力系数存在几个数量级的差异［3］，而双流体模型对微元网格内部的平均化处理却忽略了该非均匀结构的影响。

杨宁等［4］利用EMMS 模型考虑了宏观非均匀结构对计算网格的曳力修正，更准确地预测了快速流化床的颗粒夹带量。

格子玻尔兹曼方法格子玻尔兹曼方法（Lattice Boltzmann Method，简称LBM）是一种基于微观粒子动力学的计算流体力学方法，它是由Lattice Gas Automata（LGA）经过演化和发展而来的。

LBM是一种离散的方法，它通过在空间网格上模拟分子碰撞和传输过程来描述流体的宏观运动。

与传统的有限差分法、有限体积法相比，LBM具有计算效率高、并行性好、适应复杂边界条件等优点，因此在流体力学领域得到了广泛的应用。

LBM的基本思想是将流体系统离散化，将连续的流体宏观运动转化为离散的微观碰撞和传输过程。

在LBM中，流体被看作是由大量微观粒子组成的，这些微观粒子在空间网格上按照一定的规则进行碰撞和传输。

通过对微观粒子的运动状态进行统计，可以得到流体的宏观性质，如密度、速度等。

LBM的核心是格子玻尔兹曼方程（Lattice Boltzmann Equation，简称LBE），它描述了微观粒子在空间网格上的运动规律。

在LBM中，流体的宏观性质由分布函数来描述，分布函数是表示在某一时刻某一空间点上流体微观粒子的分布情况。

在每个时间步内，分布函数按照一定的规则进行碰撞和传输，通过迭代计算可以得到流体在空间网格上的演化过程。

LBM的计算过程可以并行化，因此在计算效率上具有明显的优势。

LBM的另一个优点是它对复杂边界条件的处理能力强。

由于LBM是基于离散网格的方法，因此可以比较容易地处理复杂的边界条件，如曲面边界、移动边界等。

这使得LBM在模拟复杂流体系统时具有一定的优势。

除此之外，LBM还有一些其他的优点，如对多相流、多孔介质流动等复杂流体现象的模拟能力强，对于非稳态流动和湍流流动的模拟也有一定的优势。

总之，格子玻尔兹曼方法作为一种新兴的计算流体力学方法，具有诸多优点，逐渐得到了流体力学领域的广泛关注和应用。

随着计算机硬件性能的不断提升，LBM的应用前景将更加广阔，相信它会在流体力学领域发挥越来越重要的作用。

第46卷 第1期华北理工大学学报(自然科学版)V o l .46 N o .12024年01月J o u r n a l o fN o r t hC h i n aU n i v e r s i t y o f S c i e n c e a n dT e c h n o l o g y (N a t u r a l S c i e n c eE d i t i o n )J a n .2024收稿日期:2023-03-17 修回日期:2023-12-18基金项目:国家自然科学基金项目(32070669)㊂ 第一作者:陈梦涵,女,硕士研究生㊂研究方向:应用数理统计㊂E -m a i l :2094953965@q q .c o m. 通讯作者:王希胤,男,博士,教授㊂研究方向:生物信息学㊂E -m a i l :w a n g x i y i n @v i p.s i n a .c o m. D O I :10.3969/j.i s s n .2095-2716.2024.01.013文章编号:2095-2716(2024)01-0103-08K d V 方程的格子B o l t z m a n n 模型求解陈梦涵,王希胤,李金(华北理工大学理学院,河北唐山063210)关键词:K d V 方程;D 1Q 5模型;格子B o l t z m a n n 方法摘 要:浅水波模型被广泛地用于模拟水波传播的动力学行为㊂很多问题,如强非线性问题㊁非平衡问题㊁实际应用中发生的问题等,使得传统的理论研究手段通常无能为力㊂文章首先给出了格子B o l t z m a n n 方法(L B M )的基本理论,然后利用经典的一维五速度(D 1Q 5)的离散速度模型,给出K o r t e w e g -d eV r i e s (K d V )方程中含有修正项的格子B o l t z m a n n (L B )模型推导公式,最后进行数值模拟,将K d V 方程的精确解和模拟解进行比较,然后验证修正模型的精确性㊂实验结果表明,用格子B o l t z m a n n 方法对K d V 方程进行求解,其模拟解和精确解吻合度较高㊂中图分类号:O 212.1 文献标识码:A波动现象是自然界最普遍的现象之一,对于水波的研究也一直是科学和工程研究领域的重要课题㊂尽管波动问题所涉及的领域不同,但是描述波动现象的方程却是相同的㊂由于水波千姿百态,用肉眼就可以观察到,因此很早就引起了人们的注意,可以说是人们最为熟悉的一种波㊂波动是物质运动的重要形式,广泛存在于自然界㊂波动中被传递的物理量的扰动和振动有多种形式,例如,弦线中的波㊁空气或固体中的声波㊁水波㊁电磁波,等等㊂为了更加具体地研究各种波动,就产生了各种形式的波动方程,因此,浅水波方程也成为重要的研究对象之一㊂而K o r t e w e g -d eV r i e s (K d V )方程是1895年由荷兰数学家科特韦格(K o r t e w e g )和德弗里斯(d eV r i e s )在研究浅水中小振幅长波运动时共同发现的一种单向运动浅水波偏微分方程㊂在求解偏微分方程地过程中,我们经常用到的数值计算方法有:有限元法(F i n i t eE l e m e n tM e t h o d s),有限差分法(F i n i t eD i f e r e n c eM e t h o d s ),有限体积法(F i n i t eV o l u m e M e t h o d s )和格子B o l z m a n n (L B M )方法等㊂其中,有限差分法虽然相对其他三种方法而言简便易行,而且有丰富多样的离散方法,但是它在求解问题时对求解区域的适应性比较差㊂有限元法虽然采用的网格剖分更加灵活,从一定程度上讲对求解区域具有更强的适应性,但是它在求解间断问题时会受到很大的限制,达不到有限差分法的效果㊂而有限体积法可以被视为是上述两种方法的结合,虽然能够充分利用有限元网格灵活性和克服差分法对网格适应性差的缺陷,但是数值实验较难进行[1]㊂作为一种新兴的数值模拟方法,L B M 基于B o l t z m a n n 方程的离散,是一种自下而上的求解方法㊂它描述了微观粒子的碰撞和迁移,利用分布函数(一种概率密度分布函数)来确定粒子的分布,即分布函数描述了流体的宏观运动㊂近年来,由于L B M 具有计算简便㊁良好的并行性㊁处理不规则的复杂边界容易且对于源项的考虑简单等诸多优势,已经自然而然地发展成为了求解浅水波方程的一种新方法[2]㊂在以往的研究中,英国利物浦大学的教授Z h o u [3]较为全面地阐述了浅水波方程的L B M 理论,包括外力的不同处理格式㊁湍流模型的构造㊁多种边界条件的处理方法以及对于许多经典浅水波问题的验证㊂中山大学环境科学与工程学院的L i 和H u a n g [4]进行了对流-扩散方程与浅水波方程耦合的研究,并采用L B 的多松弛模型和双松弛模型分别对流场和污染物场进行了模拟㊂文献[5]提出了一类粘性浅水方程的晶格B o l t z m a n n (L B )模型,该模型采用源项的二阶矩来恢复控制方程中的粘性,并消除C h a p m a n -E n s k o g 分析过程中产生的附加误差㊂文献[6]建立了一种适用于浅水方程的晶格玻尔兹曼模型(L A B S W E ),它用源项如床面坡度,床面摩擦力来求解方程㊂通过求解定常和非定常流动问题,验证了该模型的有效性㊂鉴于以上背景,文章首先给出了格子B o l t z m a n n 方法(L B M )的基本理论,然后利用经典的一维五速度(D 1Q 5)的离散速度模型,给出K d V 方程中含有修正项的格子B o l t z m a n n (L B )模型推导公式,最后进行数值模拟,将K d V 方程的精确解和模拟解进行比较,然后验证修正模型的精确性㊂实验结果表明,用格子B o l t z m a n n 方法对K d V 方程进行求解,其模拟解和精确解吻合度较高㊂1方法格子玻尔兹曼方法(L a t t i c eB o l t z m a n n M e t h o d)是一种基于微观介观的流体力学计算方法,适用于二维或三维流体流动问题的模拟[7]㊂图1是格子玻尔兹曼方法的流程图:其主要思想是离散化流体的分布函数,通过对分布函数的演化来模拟流体的运动㊂近年来,L B M 由于计算简单㊁并行性好㊁易于处理复杂不规则的边界及能简单方便地考虑源项等优势,已经发展成为求解浅水波方程(S W E s )的一种新方法㊂下面是格子玻尔兹曼方法的求解步骤:图1 L B M 方法流程图(1)确定格子和速度模型:首先需要确定流场的离散化格点和速度模型㊂通常情况下,将流体分成若干个小区域,每个小区域都对应一个格点,格点上有一组离散的速度向量㊂(2)定义分布函数:为了描述格点上流体的状态,引入一个分布函数g ,用来表示在每个格点上,每个速度方向上的粒子数密度㊂它是时间和位置的函数,通常用离散的速度和离散的时间步长表示㊂(3)离散B o l t z m a n n 方程:基于B o l t z m a n n 方程,对分布函数进行离散化,得到离散化的B o l t z m a n n 方程,它描述了分布函数的演化过程㊂在格子玻尔兹曼方法中,B o l t z m a n n 方程可以看成是一个简单的微分方程,其左侧是分布函数的时间导数,右侧是一个碰撞项和一个弛豫项,用于描述粒子之间的相互作用和粒子与流体之间的相互作用㊂401 华北理工大学学报(自然科学版) 第46卷(4)离散碰撞项和弛豫项:将碰撞项和弛豫项进行离散化,得到离散化的碰撞算子和弛豫算子㊂碰撞算子用于描述粒子之间的相互作用,而弛豫算子用于描述粒子与流体之间的相互作用㊂(5)迭代求解:通过迭代求解离散化的B o l t z m a n n 方程,计算出每个格点上的分布函数,从而得到流场的速度场和密度场㊂(6)计算宏观量:根据格点上的分布函数,可以计算出宏观量,如速度㊁密度㊁压力等㊂(7)处理边界条件:对于边界处的格点,需要根据具体的物理问题设置边界条件㊂(8)模拟结束:当达到预设的模拟时间或达到收敛条件时,模拟结束㊂2模型简介2.1 离散速度模型L B M 中一维离散速度模型最常见的是D 1Q 3模型和D 1Q 5模型[8],具体如下:(1)对于D 1Q 3模型(见图2),模型参数如下: c ң=c 01-1[], c s =c3, ωi =2/3,c i2=01/6,c i2=c 2{(1)图2 D 1Q 3离散速度模型其中,ωi 为权重系数,c =Δx /Δt 为粒子迁移速度,c s 是与当地声速相关的量㊂(2)对于D 1Q 5模型(见图3),模型参数如下: c ң=c 01-12-2[],c s =c (2) ω0=12,ω1=ω2=16,ω3=ω4=112(3)图3 D 1Q 5离散速度模型其中,ωi 为权重系数,c =Δx /Δt 为粒子迁移速度,c s 是与当地声速相关的量㊂2.2 K d V 方程将L B 模型应用于K d V 方程中,需要将K d V 方程离散化成网格上的方程组,然后通过L B 模型求解这个方程组㊂具体来说,L B 模型中的速度分布函数被定义为格点上的波高,通过计算速度分布函数在不同时间和空间的演化来模拟K d V 方程的行为㊂与传统的有限差分法和有限元法相比,L B 模型具有计算效率高㊁适合并行计算等优点,因此在模拟非线性波等问题时得到了广泛应用㊂考虑非线性偏微分方程一般形式[9]: ∂u ∂t +αu ∂u ∂x +βu n ∂2u ∂x 2-γ∂2u ∂x 2+δ∂3u ∂x3=0(4)其中,u =u x ,t ()是物质在空间x 处和时刻t 时的密度,α,β,γ,δ为参数㊂当β=0,γ=0时,方程(4)化为K d V 方程:∂u ∂t +αu ∂u ∂x +δ∂3u∂x3=0(5)501 第1期 陈梦涵,等:K d V 方程的格子B o l t z m a n n 模型求解3模型推导采用D 1Q 5模型给出K d V 方程含有修正项的L B 模型推导,给出的演化方程为: g i x ң+c ңi Δt ,t +Δt ()-g i x ң,t ()=-1τg i x ң,t ()-g e q i x ң,t ()()+Δt h i x ң,t ()(6)其中,τ表示松弛时间,c ңi 为沿着流动方向的单位速度矢量,gi x ң,t ()表示在t 时刻,位于x ң处沿着离散速度方向c ңi 运动的粒子分布函数,Δt h i x ң,t ()为修正项㊂根据参考文献[10-12],将修正函数h i x ң,t ()和平衡态分布函数g e qi x ң,t ()分别定义为: h i x ң,t ()=λ1,i u 2+λ2,iu n +1,i =0~4(7) g e qix ң,t ()=ρ1u (8)其中,λ1,i ,λ2,i 和p i 为调整参数㊂宏观变量u x ,t ()定义为[13-15]: u =ðigi (9)为了得到稳定的宏观变量u ,假设分布函数g i 也处于平衡状态,且有: u =ðig e qi =ði g 0()i =ðigi (10)由(10)可得, ði gn ()i=0,n ⩾1(11)为了能够恢复到宏观方程,平衡态分布函数g e q i 和修正函数需要满足下面几个约束条件:ðic i g e qi()=0,ðic 2i g e q i ()=c 2λu ,ðic3i g e qi ()=c 3ηu (12) ði hi=0,ði c i h i ()=c λ1u 2+c λ2u n +1,ðic2i h i ()=0(13)使用多尺度分析将方程恢复到宏观方程,引入1个离散的时间尺度和3个连续的时间尺度,其具体表达形式为:t 0,t 1=εt ,t 2=ε2t ,t 3=ε3t (14)对分布函数g i 和时间导数进行Ch a p m a n -E n s k o g 展开,可得: g i =g 0()i +εg1()i +ε2g 2()i +ε3g 3()i +O ε4()(15)∂∂t =∂∂t 0+ε∂∂t 1+ε2∂∂t 2+ε3∂∂t 3+O ε4()(16)其中,ε表示任意小的参数,在宏观方程的推导过程中,不妨假设ε=Δt ,将(6)式的左边对时间和空间进行泰勒展开,并保留Οε4()项,可得 ε∂g i ∂t +ε∂∂x c i g i ()+12ε2∂2g i ∂t 2+ε2∂2∂t ∂x c i g i ()+12ε2∂2∂x 2c 2i g i ()+16ε3∂∂t +c i ∂∂x æèçöø÷3g i =-1τg i -g e qi ()+εh i +Οε4()(17)将(15)和(16)代入(17)式中得,并对比左右两边可得ε的同阶项:可以得出,O ε()系数: g e qi =g0()i (18)O ε1()系数:601 华北理工大学学报(自然科学版) 第46卷∂∂x c i g e qi ()=-1τg 1()i +h i (19)O ε2()系数: ∂g e q i ∂t 1+∂∂x c i g 1()i ()+12∂2∂x 2c 2i g e qi ()=-1τg 2()i (20)O ε3()系数: ∂g e q i ∂t 2+∂g 1()i ∂t 1+∂∂x c i g 2()i ()+∂2∂x ∂t 1c i g e qi ()+12∂2∂x 2c 2i g 1()i ()+16∂3∂x 3c 3i g 1()i ()=-1τg 3()i (21)结合约束条件(12)和(13),将方程(19)两边分别乘以c i 和c 2i 后并对i 求和得:ði c i g 1()i =τði c i h i -τ∂∂x ðic 2i ge q ()i ()=c τλ1u 2+λ2u n +1()-c 2τλ∂u ∂x (22) ði c 2i g 1()i =τðic 2i h i -τ∂∂x ði c 3i g e q ()i ()=-c 3τη∂u ∂x (23)同理,结合约束条件(12)和(13),将方程(20)式两边分别乘以c i 后并对i 求和,得出:ðic i g 2()i=τ∂∂t 1ðic 2i g e q i ()+∂∂x ði c 2i g 1()i ()+12∂2∂x 2ði c 3i g e q ()i ()æèçöø÷ =c 3τ2η∂2u ∂x 2-12c 3τη∂2u ∂x 2=c 3ητ2-12τæèçöø÷∂2u ∂x2(24)结合(7)(8)(9)和(19),将方程(16)两边对i 求和,得出: ∂u ∂t 1+2c τλ1u ∂u ∂x +n +1()c τλ2u n ∂u ∂x +c 2λ12-τæèçöø÷∂2u ∂x 2=0(25)同理,结合(10)(11)(12)和(22)~(24),将方程(20)两边对i 求和,得出: ∂u ∂t 2+c 3ητ2-τ-16æèçöø÷∂3u ∂x 3=0(26)将3.19()ˑε+3.20()ˑε2,可得: ∂u ∂t 1+2c τλ1εu ∂u ∂x +n +1()c τλ2εu n ∂u ∂x +c 2λε12-τæèçöø÷∂2u ∂x 2+c 3ηε2τ2-τ-16æèçöø÷∂3u ∂x 3=0(27)将方程(27)和(4)对比可得: α=2c τλ1ε,β=n +1()c τλ2ε(28) γ=c 2λ12-τæèçöø÷ε,δ=c 3ητ2-τ-16æèçöø÷ε2(29)其中,τ=12+112+δε2ηc 3,λ=γεc 2τ-1/2(),λ1=α2τεc ,λ2=βn +1()τεc 可以得出,方程(27)就是一维K D V 方程的L B 模型㊂由方程(7)和(13)式,得出修正函数h i 为: h 0=-12λ1u 2-12λ2u n +1h 1=h 2=13λ1u 2+13λ2u n +1h 3=16λ1u 2+16λ2u n +1h 4=-13λ1u 2-13λ2u n +1ìîíïïïïïïïïïï(30)701 第1期 陈梦涵,等:K d V 方程的格子B o l t z m a n n 模型求解联立方程(10)和(12)式,可得平衡态分布函数f e q i 为:g e q 0=η2+6ρ4-λ+1æèçöø÷u =ρ0u g e q 1=λ-η2-4ρ4æèçöø÷u =ρ1u g e q2=λ2-η6-4ρ4æèçöø÷u =ρ2u g e q 3=η6+ρ4æèçöø÷u =ρ3u g e q 4=ρ4u ìîíïïïïïïïïïïïï(31)4数值模拟将L B 模型应用于K d V 方程中,需要将K d V 方程离散化成网格上的方程组,然后通过L B 模型求解这个方程组㊂具体来说,L B 模型中的速度分布函数被定义为格点上的波高,通过计算速度分布函数在不同时间和空间的演化来模拟K d V 方程的行为㊂与传统的有限差分法和有限元法相比,L B 模型具有计算效率高㊁适合并行计算等优点,因此在模拟非线性波等问题时得到了广泛应用㊂(1)考虑如下的K d V 方程:∂u ∂t +6u ∂u ∂x +∂3u∂x3=0设置参数如下:Δx =0.001,Δt =0.00001,边界条件u 0,t ()=u 255128π,t æèçöø÷=0,t >0,初始条件u x ,0()=s i n x ,x ɪ0,255128πéëêêùûúú,该方程的精确解为: u x ,t ()=c 212s e c h 2c 212x -c 21t ()+l n c 2-l n c 3éëêêùûúú取c 1=2c ,c 2=c 3=1,得出该方程的精确解为: u x ,t ()=2c 2s e c h 2c x -4c 2t ()[],x ɪ0,255128πéëêêùûúú图4 t =0.01,L B 模拟解和精确解对比 图5 t =0.25,L B 模拟解和精确解对比模拟结果如图4㊁图5所示㊂图4和图5分别给出了t =0.01和t =0.25时刻的L B 模拟解和解析解的对比图,从图中可以看出:在t =0.25之前,模拟解和解析解吻合的程度较高,但是随之时间的推移,模拟解与解析解存在一定的偏离,这主要是原因有扰动项O ε4(),它在一定程度会对孤子高度㊁速度以及形状有影801 华北理工大学学报(自然科学版) 第46卷响,且当t >0.25时,方程的模拟解和精确解差别较大㊂(2)考虑如下的K d V 方程:∂u ∂t +αu ∂u ∂x +δ∂3u∂x3=0其中,u =u x ,t (),u 为波动地振幅,x 为波横向传播的位移,t 为时间㊂初始条件为:u x ,0()=3A s e c h 2B x +C (),x ɪ0,2[]边界条件为:u 0,t ()=u 2,t ()=0,t >0解析解为:u x ,t ()=3A s e c h 2B x -D t +C (),x ɪ0,2[]其中,B =12αAδ,D =αA B .设置参数如下:Δx =0.001,Δt =5ˑ10-4,α=1,δ=4.84ˑ10-4,模拟结果如图(6)和图(7)所示. 图6 t =0.0005时,模拟解与解析解对比 图7 t =0.0020时,模拟解与解析解对比图6和图7分别给出了t =0.0005和t =0.002时刻的L B 模拟解和解析解的对比图,通过对该方程的模拟结果进行分析,发现在t =0.002之前,模拟解和解析解吻合的程度较高,但是随着时间的推移,模拟解与解析解存在一定的偏离,主要的原因是含有扰动项O ε4(),当然也可能是该类波的传播速度极快,长时间模拟就会产生偏差㊂表1给出了该方程在不同时刻的误差㊂表1 方程在不同模拟时刻的误差比较时刻/tt =0.0005t =0.0010t =0.0015t =0.0020L ɕ3.2371ˑ10-43.8721ˑ10-42.0518ˑ10-33.6037ˑ10-3L 27.4783ˑ10-55.8732ˑ10-56.4976ˑ10-47.1335ˑ10-4G R E3.8516ˑ10-44.7039ˑ10-42.5450ˑ10-34.6530ˑ10-3 从表1可以看出,在L B 模型下得到的模拟解和解析解非常逼近,无论是L ɕ误差,还是均方根误差L 2和整体相对误差G R E ,其两者之间的误差数量级都达到了,说明该数值结果是比较理想的㊂5结论现在,微分方程无处不在,各个科学领域的研究都伴随着微分方程模型㊂由于实际生活中的微分方程模型形式日趋复杂,为了与实际问题相匹配,微分方程解的形式越来越多样化㊂本文对两个特殊的K D V 方程,利用格子B o l t z m a n n 模型求解并与其精确解进行比较,得出使用格子B o l t z m a n n 方法对非线性偏微分方程求解取得了较好的效果㊂在未来的工作中,将尝试继续改进格子B o l t z m a n n 模型,并对更加复杂的偏微分901 第1期 陈梦涵,等:K d V 方程的格子B o l t z m a n n 模型求解011华北理工大学学报(自然科学版)第46卷方程或者浅水波方程进行模拟㊂希望本文可以为其他学者在求解偏微分方程方面的研究工作提供一定的参考价值㊂参考文献:[1]张海军.求解浅水波方程的熵稳定格式研究[D];西安:长安大学,2018.[2]陈文文,张文欢,汪一航,等.浅水波方程的一类改进的格子B o l t z m a n n模型[J].宁波大学学报(理工版),2020,33(01):72-79.[3] R I C K,S A L MO N.T h e l a t t i c eB o l t z m a n nm e t h o d a s ab a s i s f o r o c e a n c i r c u l a t i o nm o d e l i n g[J].J o u r n a l o fM a r i n eR e s e a r c h,1999,57(3).[4]冯士德,赵颖,茑原道久,等.旋转流场中的格子波耳兹曼模型[J].地球物理学报,2002,45(2):170-175.[5] Y U L,Z C AC,X G D,e t a l.A l a t t i c eB o l t z m a n nm o d e l f o r t h e v i s c o u s s h a l l o ww a t e r e q u a t i o n sw i t h s o u r c e t e r m s[J].J o u r n a l o fH y-d r o l o g y,2021.[6] Z H O UJG.L a t t i c eB o l t z m a n nm o d e l f o r t h e s h a l l o w w a t e r e q u a t i o n s[J].C o m p u t e rM e t h o d s i nA p p l i e d M e c h a n i c s a n dE n g i n e e r i n g,2002,191(32):3527-39.[7]张宗宁.基于格子B o l t z m a n n方法求解若干非线性偏微分方程[D];银川:北方民族大学,2022.[8]戴厚平,郑洲顺,段丹丹.一类偏微分方程的格子B o l t z m a n n模型[J].计算机工程与应用,2016,52(3):21-26.[9]王慧敏.非线性偏微分方程中孤波解的格子B o l t z m a n n模拟[D];长春:吉林大学,2014.[10]何郁波,林晓艳,董晓亮.应用格子B o l t z m a n n模型模拟一类二维偏微分方程[J].物理学报,2013,62(19):290-296.[11]乐励华,高云,刘唐伟.偏微分方程求解的一种新颖方法--格子B o l t z m a n n模型[J].大学数学,2011,27(03):75-82.[12]张春泽,程永光,李勇昌.二维浅水波方程格子B o l t z m a n n算法的G P U并行实现[J].水动力学研究与进展A辑,2011,26(2):194-200.[13]赫万恒,钱跃竑.浅水波方程的格子B o l t z m a n n模拟[C].中国力学学会北方七省市区第十三届学术大会论文集.郑州.2010:42-45.[14]赖惠林,马昌凤.非线性偏微分方程的高阶格子B G K模型[J].中国科学(G辑:物理学力学天文学),2009,39(07):913-922.[15]施卫平,胡守信,阎广武.用格子B o l t z m a n n方程模拟浅水波问题[J].力学学报,1997,(05):7-11.S o l u t i o n M e t h o d f o rK D VE q u a t i o nb y L a t t i c eB o l t z m a n n M o d e lC H E N M e n g-h a n,WA N G X i-y i n,L I J i n(C o l l e g e o f S c i e n c e,N o r t hC h i n aU n i v e r s i t y o f S c i e n c e a n dT e c h n o l o g y,T a n g s h a nH e b e i063210,C h i n a)K e y w o r d s:K D Ve q u a t i o n;D1Q5m o d e l;l a t t i c eB o l t z m a n nm e t h o dA b s t r a c t:S h a l l o w w a t e rw a v em o d e l i sw i d e l y u s e d t o s i m u l a t e t h e d y n a m i c b e h a v i o r o fw a t e rw a v e p r o p a-g a t i o n i n o c e a n a n d a t m o s p h e r e f i e l d.M a n y p r o b l e m s,s u c h a s s t r o n g n o n l i n e a r p r o b l e m s,n o n-e q u i l i b r i u m p r o b l e m s,p r o b l e m s i n p r a c t i c a l a p p l i c a t i o n s,m a k e t h e t r a d i t i o n a l t h e o r e t i c a l r e s e a r c hm e t h o d s a r e u s u a l l y p o w e r l e s s.I n t h i s p a p e r,t h e b a s i c t h e o r y o f l a t t i c eB o l t z m a n nm e t h o d(L B M)w a s f i r s t l y g i v e n.T h e n,t h ef o r m u l a o f l a t t i c eB o l t z m a n n(L B)m o d e lw i t hc o r r e c t i o n t e r mi nK o r t e w e g-d eV r i e s(K d V)e q u a t i o nw a sg i v e nb y u s i n g t h e c l a s s i c a l o n e-d i m e n s i o n a l f i v e-v e l o c i t y(D1Q5)d i s c r e t ev e l o c i t y m o d e l.F i n a l l y,t h e a c-c u r a t e s o l u t i o no f t h eK d Ve q u a t i o nw a s c o m p a r e dw i t ht h e s i m u l a t e ds o l u t i o n,a n d t h e nt h ea c c u r a c y o f t h em o d i f i e dm o d e l i s v e r i f i e d.T h e e x p e r i m e n t a l r e s u l t s s h o wt h a t t h e l a t t i c eB o l t z m a n nm e t h o d i s u s e d t o s o l v e t h eK d Ve q u a t i o n,a n d t h e s i m u l a t i o n s o l u t i o n i s i n g o o d a g r e e m e n tw i t h t h e e x a c t s o l u t i o n.。

基于格子波尔兹曼方法的回热器数值模拟夏宇栋;陈曦;马诗旻;张华【摘要】Based on Lattice Boltzmann Method ( LBM ) , the flow field in the microstructure of mesh screen and etched foil regenerators were simulated. The velocity and pressure drop were obtained by LBM. The simulation results show that velocity field in etched foil regenerator is better distributed than that of mesh screen. And the etched foil regenerator has less resistance coefficient than that of mesh screen regenerator. The simulation results basically agree with the experiments.%利用格子玻尔兹曼方法,直接对蚀刻薄片和层叠丝网回热器的微观结构流场进行了模拟.得到了两种回热器填料的微观流场和两端的压差.模拟结果显示,当回热器的直径、水力直径和填充率相近情况下,不同流速下蚀刻薄片卷裹式回热器的稳态阻力系数均比层叠丝网回热器小.稳态阻力系数的模拟变化趋势与实验一致.【期刊名称】《低温工程》【年(卷),期】2012(000)005【总页数】5页(P41-45)【关键词】格子波尔兹曼方法;回热器;流阻系数;数值模拟【作者】夏宇栋;陈曦;马诗旻;张华【作者单位】上海理工大学能源与动力学院上海200093;上海理工大学能源与动力学院上海200093;上海理工大学能源与动力学院上海200093;上海理工大学能源与动力学院上海200093【正文语种】中文【中图分类】TB651回热器是低温制冷机的关键部件，也是影响制冷机性能的最重要因素。

格子Boltzmann方法三种边界格式的对比分析刘连国;杨帆;王宏光【摘要】采用格子Boltzmann方法(Lattice Boltzmann Method- LBM)对二维顶盖驱动方腔流动进行数值模拟.在计算中分别使用半步长反弹、非平衡反弹、以及非平衡外推三种边界处理格式,并得到了不同格式对应的流线分布,流函数最小值、涡心坐标、几何中心线速度分布等.通过将所得结果与基准解进行比较,就三种边界格式的计算效率,计算精度、以及计算稳定性等方面进行了讨论和分析,为LBM计算中边界格式的选择提供了有益的参考.【期刊名称】《机械研究与应用》【年(卷),期】2012(000)001【总页数】5页(P18-22)【关键词】格子Boltzmann方法;边界处理格式;半步长反弹格式;非平衡反弹格式;非平衡外推格式【作者】刘连国;杨帆;王宏光【作者单位】上海理工大学能源与动力工程学院,上海200093;上海理工大学能源与动力工程学院,上海200093;上海理工大学能源与动力工程学院,上海200093【正文语种】中文【中图分类】O357.11 引言格子Boltzmann方法(LBM)是近年来迅速发展的一种新型数值计算方法。

边界条件的处理是LBM实施中一项非常关键的内容。

实际计算表明:选取不同的边界条件会对数值计算的精度、稳定性以及效率产生很大影响。

作为LBM的一个基本问题，边界条件的处理一直是流体力学一个重要的研究方面。

根据边界条件的类型，可将之分为两类:压力边界和速度边界［1］，其中的速度边界又可细分为:平直边界和曲面边界。

笔者从经典的流体力学问题二维顶盖方腔流模拟入手，对三种平直边界格式进行对比和分析，为LBM计算中边界格式的选择提供了有益的参考。

2 二维九点格子Boltzmann模型目前最常用的格子Boltzmann模型为LBGK模型，通过引入“单一弛豫时间”来简化Boltzmann方程中碰撞项的计算［2］。

九点格子LBGK模型的演化方程为:式中:(x，t)是在t时刻、x处的平衡态分布函数;τ为单一弛豫时间因子;eα为网格点各方向上的粒子速度。

格子玻尔兹曼算法
格子玻尔兹曼算法是一种基于微观粒子运动的计算流体力学方法，它可以用来模拟流体的运动和传输过程。

该算法的核心思想是将流体分成许多小的格子，然后在每个格子内模拟流体粒子的运动和相互作用，从而得到整个流体的宏观运动状态。

格子玻尔兹曼算法的基本原理是通过离散化的方式来模拟流体的微观运动。

在每个格子内，流体粒子的运动状态可以用一个分布函数来描述，该函数包含了流体粒子在不同速度下的密度和速度信息。

通过对分布函数的离散化和更新，可以得到流体的宏观运动状态，如速度、密度和压力等。

格子玻尔兹曼算法的优点是可以处理复杂的流体运动和传输过程，如湍流、多相流和热传导等。

同时，该算法具有高效、可扩展和易于并行化等特点，可以在大规模计算机集群上进行高性能计算。

然而，格子玻尔兹曼算法也存在一些挑战和限制。

首先，该算法需要对流体的微观运动进行离散化，因此需要选择合适的离散化方法和参数，以保证模拟结果的准确性和稳定性。

其次，该算法需要进行大量的计算和存储，因此需要高性能计算机和存储系统的支持。

最后，该算法在处理复杂流体问题时，需要考虑多种物理过程的相互作用，因此需要进行多物理场的耦合和协同计算。

格子玻尔兹曼算法是一种重要的计算流体力学方法，它可以用来模
拟各种复杂的流体运动和传输过程。

随着计算机技术的不断发展和进步，该算法将在更广泛的领域得到应用和发展。

传热学格子玻尔兹曼方法计算方法的特点摘要本文讨论了传热学中的格子玻尔兹曼方法，并分析了这一计算方法的特点。

首先，我们介绍了传热学的基本概念和研究背景。

然后，我们详细解释了格子玻尔兹曼方法的原理和模拟过程。

接着，我们探讨了该方法的特点，包括计算效率、模拟精度和适用范围等。

最后，我们总结了格子玻尔兹曼方法在传热学中的应用前景，并提出了进一步研究的方向。

1.引言传热学是研究能量从一个物体传递到另一个物体的学科。

在工程领域中，传热问题经常出现在热流体系统的设计和优化中。

传热过程涉及热传导、对流和辐射等多种传热机制，准确模拟传热过程对于工程实践和科学研究具有重要意义。

格子玻尔兹曼方法（L a tt ic eB ol tz ma nnM e th od,L BM）是一种基于微观颗粒模拟传输过程的计算方法，近年来在传热学领域得到了广泛应用。

与传统的求解传热方程的数值方法相比，格子玻尔兹曼方法通过模拟颗粒在格子上的运动来描述流体的宏观行为，具有更高的计算效率和更灵活的模拟能力。

2.格子玻尔兹曼方法原理格子玻尔兹曼方法基于玻尔兹曼方程和格子自动机理论，通过在一个规则的网格上模拟微观颗粒的运动来模拟流体的运动。

格子玻尔兹曼方法的基本原理是将流体分割成一系列小的正方体，每个正方体称为格子。

在每个格子中，通过对流、碰撞和反弹等过程来模拟颗粒之间的相互作用。

格子玻尔兹曼方法的模拟过程可以分为以下几个步骤：1.确定模拟区域的网格分布和流体的边界条件。

2.初始化流体的宏观和微观状态，在格子中随机分布将流体颗粒的速度和密度初始化为一定状态。

3.对于每个时间步长，根据碰撞和对流过程更新格子中流体颗粒的状态。

4.根据流体颗粒的状态计算宏观流体变量，如流速和压力等。

5.重复步骤3和4，直到达到设定的模拟时间。

3.格子玻尔兹曼方法特点格子玻尔兹曼方法具有以下几个特点：3.1计算效率高格子玻尔兹曼方法在模拟复杂流体系统时具有较高的计算效率。

格子boltzmann方法的理论及应用
格子波尔兹曼方法(Grid Boltzmann Method, GBM)是一种非离散化处理方法，其基本
思想是在空间上采用格点，并建立格点微分方程组来解决复杂流体或者其他相关物理问题. GBM以较少的计算量就可达到快速、精确求解流体动力学问题，而且将空间和时间分离，
大大减少计算量和存储量，可以说是比传统有限元技术和有限差分技术更加有效的一种方法.
格子波尔兹曼方法的具体原理是：格子波尔兹曼方法是将空间上的解释解划分成一系
列的蒙特卡洛格子点，这样可以以非离散化处理。

针对与流体物理仿真相关的变量，以格
点位置为基底，可以使用波尔兹曼分布Y(v)来描述，将原本复杂的多体相互作用模型转化为简单的蒙特卡洛定值模型，由此通过空间离散的方式可以求解波尔兹曼方程；具体的应
用也很广泛，可以应用在流体动力学中，可用来模拟很多液体问题，比如湍流传播和燃烧
等方面；在地形风化中可以用来模拟流域洪水演变和地形演化、土壤流失问题；在水质污
染领域，可以用来模拟河流污染物质运行规律；在非牛顿流体中，可用来模拟非牛顿流体
动力学问题；在金属粒子、微粒或者多组分液体中，可用来模拟粒子间相互作用，甚至可
以应用在非弹性波中进行数值模拟.
格子波尔兹曼方法因其独特的优越性深受广泛重视，在国内外都有大量的研究，结合
其他的数值方法，用于模拟复杂的流体物理系统，改善计算效率，提高建模的准确性。

GBM具有更快的计算速度和精度优势，在现代的科学技术领域有着广泛的应用，如流体动
力学，地形风化，水质污染等问题。

该方法不仅可用作模拟计算复杂流体运动，而且可以
用于半定常及强力学分析中。

有限体积格子Boltzmann方法的算法改进及性能分析武频;曹啸鹏;尚伟烈;郑德群;高升【期刊名称】《计算机应用研究》【年(卷),期】2012(29)10【摘要】有限体积格子Boltzmann方法(LBM)能够将标准LBM的应用范围扩展到非结构网格,但是比起标准的LBM这个方法需要更多的内存用量和计算量.针对此问题采用了优化计算顺序、简化计算方程的方法对有限体积LBM算法进行改进.科学的分析和实验的结果表明,改进后的算法能够在不增加计算量的基础上减少内存用量,在一些情况下还可以大量减少计算时间.%Finite volume lattice Boltzmann method ( LBM) can extend the standard LBM to unstructured mesh, but compared with standard LBM this method suffers from higher memory consumption and poorer computational performance. In order to solve this problem, the improvement process used the methods of optimizing evaluation order and simplifying calculating equation. Scientific analysis and experimental results demonstrate that the improved algorithm results in lower memory usage without additional computation, and in some conditions it reduces much computation.【总页数】4页(P3706-3709)【作者】武频;曹啸鹏;尚伟烈;郑德群;高升【作者单位】上海大学计算机工程与科学学院,上海200072;上海大学计算机工程与科学学院,上海200072;上海大学计算机工程与科学学院,上海200072;上海大学计算机工程与科学学院,上海200072;上海大学计算机工程与科学学院,上海200072【正文语种】中文【中图分类】TP391【相关文献】1.异构平台下格子Boltzmann方法实现及性能分析 [J], 张丹丹;徐莹;徐磊2.基于格子Boltzmann方法和有限体积法的方柱绕流特性对比分析 [J], 史冬岩;李红群;王志凯3.基于非结构化网格的高可扩展并行有限体积格子Boltzmann方法 [J], 徐磊;陈荣亮;蔡小川4.有限体积格子Boltzmann方法用于近空间连续流区绕流模拟 [J], 皮兴才;李志辉;彭傲平;张子彬5.格子Boltzmann并行程序的优化与性能分析 [J], 赵鹏;张丹丹;汪鲁兵;田振夫;钱跃竑因版权原因，仅展示原文概要，查看原文内容请购买。

格子boltzmann法
格子波尔兹曼法（Grid-Based Boltzmann Method）是用于计算复杂系统的一种数值模拟方法，该方法基于玻尔兹曼方程，采用格子划分的非总熵方案计算分布函数所描述的传播动力学系统的平衡性质。

格子波尔兹曼方法由三个部分组成，分别是分子动力学基础、格子化方案以及格点迭代方案。

在空间上，格子波尔兹曼方法采用密度聚类格子，由于每个格子内节点之间的影响，允许改变每个节点状态。

在时间上，格子波尔兹曼方法通过欧拉法和龙格-库塔法，将弹性系统的猝灭问题转换为一个接近平衡态的迭代问题。

最终，根据初始条件和格子化方案计算本征周期、如粒子操纵力学系统中的陷阱模式等。

热格子Boltzmann法分析及应用陈杰;钱跃竑【摘要】格子Boltzmann方法(lattice Boltzmann method,LBM)是一种基于气体动理论的介观计算方法,其物理背景清晰、边界处理简单,已成功应用于等温(或无热)流动中.简要介绍现有的几种热格子Boltzmann模型,并运用几种热格子模型求解热Couette流、方腔自然对流等典型算例,对比不同热格子模型的数值稳定性、准确性、模型的计算效率等.将两种热格子模型用于多孔介质内的流动与传热问题中,对比热格子模型在处理复杂结构时的数值特性.%Lattice Boltzmann method (LBM) is a mesoscale computational method based on the gas kinetic theory. For solving Fourier-Navier-Stokes equations, the thermal lattice model has attracted much research attention. This paper compares several thermal lattice models in terms of accuracy, stability and computational efficiency. The thermal flow in pore-scale porous is also studied using different thermal lattice models.【期刊名称】《上海大学学报（自然科学版）》【年(卷),期】2012(018)005【总页数】7页(P489-495)【关键词】格子Boltzmann方法;热格子Boltzmann方法;多孔介质【作者】陈杰;钱跃竑【作者单位】上海大学上海市应用数学和力学研究所,上海200072;上海大学上海市应用数学和力学研究所,上海200072【正文语种】中文【中图分类】O351格子Boltzmann方法(lattice Boltzmann method，LBM)是近20年发展成熟起来的一种数值计算方法.LBM基于气体动理论，通过分布函数的演化获得宏观信息.作为一种简单且能处理复杂流动问题的有效数值方法［1-2］，LBM具有良好的数值稳定性、天然的并行性、简单的边界处理等优点，自出现之日起就被广泛用于多孔介质流［3］、多相流［4］、反应扩散系统［5］等诸多领域.早期的LBM只应用于等温流动(或无热流动)的模拟，但是基于这种方法具备处理复杂问题的能力以及解决传热问题的需要，研究者一直在不断地探索研究热格子Boltzmann模型，已形成了一些经过数值验证具有模拟热流动能力的热LBM［6-10］，并应用于多孔介质流动与传热、燃烧及化学反应流、湍流等问题.本研究简述了不同热格子Boltzmann模型的基本理论，并通过数值分析对比了不同热格子Boltzmann模型的计算结果及数值特性，进而用于多孔介质流动传热问题中.1 等温LBM基本原理LBM中除时间、空间被离散之外，无限维的粒子速度空间也都被离散成有限的速度序列.在标准LBM模型中，物理空间被离散成正方形(体)格子，流体粒子在格点x上碰撞并按离散速度E=［e0，e1，…，eq－1］迁移到x+eiδt格点.fi(x，t)定义为t时刻在格点x上速度为ei的粒子密度，满足如下的格子Boltzmann方程：式中为平衡态函数，ω为松弛因子.通过简单地向平衡态不断趋近的过程代替真实的复杂碰撞，即BGK(Bhatnagar-Gross-Krook)近似，所以此模型也称为LBGK 模型.平衡态分布函数的选取是LBM的关键.DnQm系列［1］中均采用式中，cs为格子声速，Wi为不同速度粒子的权重.本研究在数值模拟中均采用D2Q9模型.宏观密度和速度分别定义为2 热格子Boltzmann模型现有的热格子Boltzmann模型通常可以分为两大类：第一类是流场温度场耦合统一求解的模型，如多速格子Boltzmann模型(multi-speed LBM，MSLBM)、熵格子Boltzmann方法(entropic LBM，ELBM);另一类则是对流场与温度场分别求解，如被动标量格子Boltzmann模型(passive scalar LBM，PSLBM)、双分布函数(double-distribution-function，DDF)模型，以及其他与传统计算流体动力学(computational fluid dynamics，CFD)结合的混合方法，如混合热格子Boltzmann方法(hybrid-thermal LBM，HTLBM).2.1 多速格子Boltzmann模型(MSLBM)多速格子Boltzmann模型是等温LBM模型的直接推广，其密度、速度、内能等均由速度分布函数的各阶速度矩得到.Qian［6］基于等温LBGK模型，提出了D1Q5，D2Q13，D3Q21，D3Q25热力学LBGK模型.在这些模型中，除了要满足等温模型的守恒条件外，还应满足能量守恒和平衡态热通量为0的条件：平衡态分布函数是Maxwell分布的截断形式：式中，Ap，Bp，Dp为待定参数，由满足的守恒条件确定.平衡态包含了速度的三阶项，离散速度也在D2Q9的基础上在主坐标轴上增加了4个速度.Qian［6］采用此模型对一维激波管、二维 Rayleigh-Benard对流进行了模拟，证明了该模型的有效性.MSLBM具有良好的物理基础，宏观方程绝对耦合，已成功模拟了一些传热现象，但只能模拟狭窄的温度范围和较小的Ma数，存在稳定性问题，限制了该模型的广泛应用.2.2 熵格子Boltzmann方法(ELBM)熵格子Boltzmann方法考虑了H定理，通过在守恒约束下最小化波尔兹曼H函数求解平衡态分布函数，由此得出的正定的分布函数保证了模型的稳定性和准确性［11］.Prasianakis等［10］将ELBM拓展到热流动问题的求解中，证实了该方法的有效性，本研究参照此方法.H函数定义为平衡态分布函数则是在满足守恒约束条件：的情况下，求H函数最小值得到的，具体形式详见文献［10］.Prasianakis等［12］采用在ELBM中加入高阶量的补偿算法，较大地提高了基于D2Q9标准格子的ELBM可模拟的温差和Ma数，但是模型实施较为复杂.2.3 双分布函数模型双分布函数模型，即存在两个分布函数：密度分布函数和内能(温度或总能)分布函数，其中密度分布函数用于模拟速度场，而内能(温度或总能)分布函数则用来模拟温度场.温度、内能或总能分布函数均通过不同的方式构造，但其演化都独立于密度分布函数.2.3.1 被动标量格子Boltzmann模型(PSLBM)被动标量格子Boltzmann模型基于如下原理：在忽略压力做的功和粘性热耗散的情况下，温度可以看作是随流体运动的一个标量，遵循对流扩散方程.由于此方程与组分浓度场的控制方程一样，于是Shan［7］提出使用两组分模型模拟单组分热流动问题：组分1模拟流体的运动;组分2模拟被动的温度场.平衡态密度函数为式中，σ表示组分，两组分共享速度，2.3.2 内能双分布函数模型内能双分布函数模型最早由He等［8］提出，其速度场仍用密度分布函数演化模拟，温度场则由内能分布函数模拟.该模型的基本思想是通过对连续Boltzmann方程进行特殊的离散得到等温LBM，如果进行同样的操作，则热LBM可以由离散内能的演化方程得到.根据内能的定义ρε=∫(ξ－u)2/2f dξ，引入内能分布函数g(r，ξ，t)=(ξ－u)2f/2，并引入新的碰撞模型，得到内能分布函数满足的演化方程：式中，q=(ξ－u)·［∂tu+(ξ·)u］.然后对演化方程离散，得到可用于数值计算的离散的分布演化方程，具体的离散过程详见文献［8］.相比于PSLBM，内能DDF的构造更具有物理基础，并包含了粘性热耗散和可压缩功.相比于MSLBM，DDF模型具有更好的数值稳定性，Pr数不受限制，因此被广泛用于各种近似不可压流体流动与传热问题.2.4 混合热格子Boltzmann模型(HTLBM)HTLBM是指使用 LBM解速度场，使用传统CFD解温度场，并通过一定的方式相互影响.这种方法利用了LBM能简单处理复杂流动问题的优势以及传统CFD在传热问题上的成熟技术，可以处理一些仅仅使用传统CFD较难解决的复杂流动传热问题.最初，Lallemand等［13］将多速多松弛模型和有限差分法(finite difference method，FDM)相结合，提出了混合模型，速度场用多松弛LBM求解，温度场采用FDM求解.本研究采用有限容积法(finite volume method，FVM)与LBM相结合的混合方法，即采用如下的FVM求解能量守恒方程：式中，S为广义源项，包括压力做的功和粘性热耗散.速度场与温度场的耦合通过在LBM中添加温度相关的外力项以及在FVM中添加广义源项S来实现.此外，普朗特数、比热容等热物性以及随温度变化的输运系数可以实现相应的调节.本研究中FVM与LBM采用同一套网格系统，FVM采用绝对稳定且具有与LBM相同精度的二阶迎风格式(second-order upwind scheme，SUS).PSLBM，DDF以及HTLBM这类模型的一个关键之处在于流场与温度场之间的耦合，其模型往往不满足气体完全状态方程，温度场对速度场的影响只是通过施加一个外力来实现.如Guo等［9］针对Boussinesq方程组，通过在密度分布函数演化方程中增加一个外力项以实现温度对流场的影响.Filippova等［14］基于HTLBM研究了小Ma数下高温燃烧，用温度场修正密度场以满足状态方程.3 计算结果及分析为了进一步对比各类模型，本研究采用ELBM，PSLBM，内能DDF模型以及HTLBM，对热Couette流、封闭方腔自然对流和多孔介质内非等温流动等问题进行了模拟对比.3.1 热Couette流模拟考虑两平板间热Couette流，上平板以速度U向右运动，下板静止，且上下平板分别保持恒温Th，Tc，且Th＞Tc.横截面温度廓线的解析形式为式中，H为平板间距离，Pr=ν/χ为普朗特数，χ为热扩散系数，Ec=U2/［Cp(Th －Tc)］为埃克特数.热Couette流中不考虑流体可压缩性的影响，而粘性耗散效应明显，因而分别运用ELBM，内能DDF模型和HTLBM对该问题进行了模拟，网格数均为64×64.模拟中Re=UH/ν=20，计算结果如图1所示.固定Pr=4，Ec分别为1，10和20的无量纲温度廓线，散点为不同方法的计算值，曲线为解析解公式(10).由图可见，三种模型都成功模拟了粘性耗散效应，且与解析解吻合得很好.本工作进一步研究了三种模型的计算效率问题.图2给出了温度残差随CPU时间的变化曲线，可见ELBM和HTLBM明显优于内能DDF模型.3.2 封闭方腔自然对流模拟封闭方腔尺寸为H(正方形边长)，左右壁面分别保持恒温Th，Tc，且Th＞Tc，上下壁面绝热，四壁面速度均为无滑移边界.方腔内充满均质空气，考虑向下的重力.描述自然对流的无量纲参数Ra数定义为图1 热Couette流温度廓线Fig.1 Temperature variation of the thermal Couette flow图2 热Couette流温度残差变化曲线Fig.2 Temperature residuals variation of the thermal Couette flow式中，β为热膨胀系数.物性满足Boussinesq假设，这里通过施加外力G=－β(T－T0)g实现温度场对速度场的影响.在方腔自然对流中，可压缩效应以及粘性耗散效应可忽略不计.从模型分析可以看出，PSLBM在这种情况下与DDF模型类似，而ELBM边界实施较为复杂.因此，本研究分别采用不包含粘性耗散效应的PSLBM和HTLBM对该问题进行了模拟，模拟中Pr=0.71，Ra数分别为104，105和106.图3和图4分别为HTLBM在不同Ra数下流动稳定后得到的流线、等温线，与以往的数值及实验结果一致.由图3可见，随着Ra数的增大，方腔中心的近似圆形的涡逐渐变成椭圆形，进而分裂成两个涡.当Ra= 106时，两个涡分别向左右壁面移动，在中心出现了第三个涡.由图4可见，随着Ra数的增大，竖直的等温线逐渐变得水平，主导的传热机理由导热变为对流.为了进一步定量考核，本研究计算了努塞尔数Nu和平均努塞尔数 Numean.表1给出了热壁面的Numean、最大Nu数Numax及相应位置的yNumax、水平中心线上最大速度vmax及相应的位置x、垂直中心线上最大速度umax以及相应的位置y.HTLBM和PSLBM求解的结果与Barakos等［15］的基准解一致.同样，本研究对HTLBM和PSLBM的计算效率进行了对比，图5所示为两种方法模拟自然方腔对流Ra=105时，速度残差随CPU时间的变化曲线.可以明显看出，两种方法中残差均呈现震荡下降趋势，且HTLBM收敛快于PSLBM，HTLBM残差收敛到10－7以下时的耗时为PSLBM的57%.图3 方腔自然对流不同Ra数的流线Fig.3 Predicted streamlines of natural convection图4 方腔自然对流不同Ra数的等温线Fig.4 Predicted temperature profiles of natural convection表1 数值解与基准解对比Table 1 Comparison of numerical results between thermal models and benchmarksRa数模型 Numean Numax(y/H) umax(y/H) vmax(x/H) PSLBM 2.247 3.538(0.141) 0.194(0.824) 0.234(0.121) Ra=104 HTLBM 2.242 3.553(0.145) 0.194(0.824) 0.234(0.121) Barakos等［16］2.2453.539(0.143) 0.193(0.818) 0.234(0.119) PSLBM4.512 7.827(0.075)0.128(0.854) 0.256(0.065) Ra=105 HTLBM 4.507 7.723(0.085) 0.134(0.854) 0.260(0.065) Barakos等［16］ 4.510 7.636(0.085) 0.132(0.859) 0.258(0.066) PSLBM 8.809 17.454(0.033) 0.079(0.852) 0.261(0.037) Ra=106 HTLBM 8.792 17.435(0.040) 0.081(0.854) 0.263(0.040) Barakos等［16］ 8.80617.442(0.037) 0.077(0.859) 0.262(0.039)图5 方腔自然对流速度残差变化曲线Fig.5 Velocity residuals variation of thenatural convection3.3 多孔介质非等温流动模拟多孔介质内部结构十分复杂，其流动传热现象也相当复杂.格子Boltzmann方法在模拟孔隙内的流体运动时可以方便地使用反弹格式处理复杂流场，因此，该方法在孔隙尺度模拟多孔介质内部复杂流动上有明显的优势及较高的计算率.对于多孔介质内流动与传热的问题，以往使用比较广泛的是PSLBM和内能DDF模型.本研究将HTLBM用于多孔介质流动与传热分析中，并与PSLBM进行了对比.本研究分析了分形多孔介质中的自然对流，分形结构采用Sierpinski地毯，依次对分形等级N=2和3的Sierpinski情况进行了模拟.无量纲控制参数Pr=0.71，Ra数分别为104，105和106，固体区域温度保持线性温度分布.图6为采用HTLBM计算N= 2分形结构内自然对流得到的流线图，图7为相应的等温线.由图可见，模拟结果与PSLBM一致，随Ra数的逐步增大，传热机理由导热主导变化为对流主导.图8为N=3，Ra=106时的流线图及等温线.由图可见，固体的增多明显地抑制了对流作用.同样对HTLBM在计算效率的问题上和PSLBM进行了对比.图9为Ra=106时两种方法模拟N=2分形结构时的速度残差曲线，此时HTLBM耗时为PSLBM的76%，仍具有优势.图6 多孔介质方腔自然对流流线(N=2)Fig.6 Predicted streamlines of porous cavity(N=2)图7 多孔介质方腔自然对流等温线(N=2)Fig.7 Predicted temperature profiles of porous cavity(N=2)图8 多孔介质方腔自然对流流线及等温线(N=3)Fig.8 Predicted streamlines and temperature profiles of porous cavity(N=3)4 结论本研究简要介绍了几种热格子Boltzmann模型(MSLBM，ELBM，PSLBM，内能DDF模型及HTLBM)，并运用不同热格子模型求解了两个典型算例以及多孔介质流动传热问题，得到如下结论.图9 多孔方腔自然速度残差变化曲线Fig.9 Velocity residuals variation of porous cavity(1)速度场温度场耦合求解的模型还需要进一步发展才能被广泛应用.(2)相比于PSLBM和DDF模型，HTLBM在保证计算精度的前提下，具有较高的计算效率.(3)数值模拟验证了HTLBM在处理多孔介质复杂结构时可行、有效，且比PSLBM 的效率高.参考文献：［1］ QIANY H，D’HUMIERESD，ttice BGK models for Navier-Stokes equation ［J］.Europhysics Letters，1992，17(6)：479-484. ［2］ QIANY H，SUCCIS，ORSZAGS A.Recent advances in lattice Boltzmann computing［M］∥ DIETRICH S.Annual reviews of computational physicsⅢ.New J ersey：World Scientific Publishing Company，1995：195-224.［3］ ZHAOC Y，DAIL N，TANGG H，et al.Numerical study of natural convection in porous media(metals) using lattice Boltzmann method (LBM) ［J］.International Journal of Heat and Fluid Flow，2010，31 (5)：925-934. ［4］严永华，石自媛，杨帆.液滴撞击液膜喷溅过程的LBM模拟［J］.上海大学学报：自然科学版，2008，14(4)：399-404.［5］李青，徐旭峰，周美莲.三维斑图形成的格子Boltzmann方法模拟［J］.上海大学学报：自然科学版，2007，13(5)：516-518.［6］ QIANY H.Simulating thermohydrodynamics with lattice BGK models ［J］.Journal of Scientific Computing，1993，8(3)：231-242.［7］ SHANX.Simulation of Rayleigh-Bénard convection using a lattice Boltzmann method［J］.Physical Review E，1997，55(3)：2780-2788. ［8］ HEX，CHENS，DOOLENG D.A novel thermal model for the latticeBoltzmann method in incompressible limit［J］.Journal of Computational Physics，1998，146 (1)：282-300.［9］ GUOZ，ZHENGC，SHIB，et al.Thermal lattice Boltzmann equationfor low Mach number flows：Decoupling model［J］.Physical Review E，2007，75 (3)：036704.［10］ PRASIANAKISN I，CHIKATAMALAS S，KARLINI V，et al.Entropic lattice Boltzmann method for simulation of thermal flows［J］.Mathematics and Computers in Simulation，2006，72(2)：179-183. ［11］ ANSUMALIS，KARLINI V，OTTINGERH C.Minimal entropic kinetic models for hydrodynamics ［J］.Europhysics Letters，2003，63(6)：798-804.［12］ PRASIANAKISN I，KARLINI ttice Boltzmann method for simulation of compressible flows on standard lattices［J］.Physical Review E，2008，78(1)：016704.［13］ LALLEMANDP，LUO L S.Theoryofthelattice Boltzmann method：Acoustic and thermal properties in two and three dimensions［J］.Physical Review E，2003，68(3)：036706.［14］ FILLIPPOVAO，HANELlD.A novellatticeBGK approach for low Mach number combustion［J］.Journal of Computational Physics，2000，158(2)：139-160.［15］ BARAKOSG，MITSOULISE，ASSIMACOPOULOSD.Natural convection flow in a square cavity revisited：Laminar and turbulent models with wall functions［J］.International Journal for Numerical Methods in Fluids，1994，18(7)：695-719.。

第41卷第1期东北电力大学学报Vol.41，H 2021年2月Journal Of Northeast Electric Power University Feb,2021D O I： 10. 19718/j.issn. 1005-2992.2021-01-0048-08求解辐射传输方程的多松弛格子-Boltzmann模型刘晓川\王存海2,黄勇、朱克勇1(1.北京航空航天大学航空科学与工程学院，北京100191,2.北京科技大学能源与环境工程学院，北京100083)摘 要：针对福射传输方程，文中提出了一种多松弛格子-B o l t z m a n n模型（m u l t i p l e-r e l a x a t i o n-t i m el a t t i c e B o ltz m a n n m o d e l).基于扩散尺度下的M a x w e l l迭代，辐射传输方程可以严格地从格子B o l t z m a n n方程推导得出•一些数值案例用来验证本文提出的多松弛（M R T)格子-B o l t z m a n n模型的数值特性.结果表明本文提出的多松弛格子-B o l t z m a n n模型可以稳定及精确地求解参与性介质中的瞬态及稳态辐射传输问题.并且，该模型具有二阶精度.关键词：辐射传输方程;格子-B o l t z m a n n方法;多松弛模型中图分类号：T K121 文献标识码：A辐射传输方程描述了辐射能量在介质中的传递，在许多科学和工程领域具有重要作用，例如大气辐 射传输[1]、光学层析成像[2]、天体物理学[3]及核工程[4]等.辐射传输方程是一个高维、复杂的积分微分 方程，辐射强度涉及波长、时间、空间和角度等,求其解析解十分困难.学者们提出发展了很多种数值方 法来求解辐射传输方程,如蒙特卡洛法[5]，离散坐标法[6]，有限体积法[7]，有限元法[8]等.近年来，利用格子-Boltzmann方法（L B M)来求解辐射传输方程吸引了许多学者的兴趣.L B M起源 于格子气自动机，已经发展成为了一种计算流体力学的有力数值工具[9].并且,L B M已经被拓展到求解 许多线性和非线性系统问题，例如声子输运[1°],波传播[11]，反应扩散,对流扩散等.相比于其他的求解辐射传输方程的数值方法，L B M不需要计算大量的光线轨迹，也不需要离散复杂的偏微分方程. L B M具有容易实现，高并行效率等优点•目前，对于利用L B M来求解辐射方程还不完善，发展完善的 L B M用于求解辐射传输方程是必要的.1^3111^等[14]假定了可调节的虚拟光速和辐射平衡条件，将L B M推广到分析参与性介质中的辐射 问题.M a等[~基于辐射流体力学，提出了一维辐射的格子-Boltzmann模型.Zhang等[16]通过采用全隐 式后项差分格式处理辐射方程中的瞬态项，将L B M扩展到求解参与性介质中的一维瞬态辐射传输. Mink等[171在将P1近似应用辐射传输方程的基础上提出了一种三维的格子-Boltzmann模型，然而此模 型仅适用于光学厚介质.Y i等[18]通过引入虚拟的扩散项，将辐射传输方程视为一种特殊的对流扩散方 程，从而提出了一种二维稳态射传输方程的格子-Boltzmann模型.W a n g等[19_将瞬态辐射传输方程处 理为双曲守恒方程，然后提出了 一■种求解瞬态辅射和中子输运的格子-Boltzm ann模型.目前，求解辐射方程的多松他的格子-Boltzmann模型还未见报道.本文提出了一种多松她格子-Bo­ltzmann 模型 （multiple-relaxation-time l a t t i c e Boltzmann mode丨）■基于扩散尺度下的 Maxwell 迭代，福射传收稿日期：2020-11-09基金项目：国家自然科学基金（s is M o o w g c te o i4)第一作者：刘晓川（1992-),男，在读博士研究生，主要研究方向：航空科学与工程通讯作者：黄勇（1974-),男，博士，教授，主要研究方向：航空科学与工程电子邮箱：liuxiaochuan@(刘晓川），wangcunhai@ustb_(王存海），huangy@(黄勇），zhukeyong@buaa.edu_cn(朱克勇）第1期 刘晓川等：求解辐射传输方程的多松弛格子-Bohmiami模型 49输方程可以严格地从格子Boltzmann方程推导得出，并且不引人任何限制和近似.本文发展的多松弛格 子-Boltzmann模型可以精确地求解参与性介质内的多维瞬态及稳态辐射传输问题.数值结果表明该模 型具有二阶精度和收敛速率.并且，相比于单松弛模型，多松弛模型具有更好的稳定性.该模型可以进一 步推广到求解参与性介质内的辐射传输问题.1福射传输方程的多松她格子-B o l t z m a n n模型1-1辐射传输方程考虑吸收、发射和散射介质内的辐射传输方程，其离散坐标形式可以写为[2°]dl(r,rr,t) +f f, v/(r；i T^)+/3(r)l(r,f r,t)=S(r,n r,t)，(1)cLdt公式中心为介质内的光速;/为辐射强度;r为位置坐标冶+屹为衰减系数;/r + V")+ 为离散方向，源项S可以表示为s(r,n r,t)^kaib(r,{T,t)+^J j i{r,i r')(p(n r\n n)w m',(2)47T m，=1公式中A为总的离散方向，=1,2,…，八^' = 1,2,…，yv;M；m'为对应方向的权重.考虑漫发射和反射壁面，边界条件可以写为I(rw,{r,t)^e wIb(rw,t)+^-^I(rw ,f T') \nw -HT'\w m' + (\ - p j r\rw J F",t) ,(3)17 <Q公式中:&为发射率;Pu，为反射率;广‘为外部人射辐射强度.1.2 多松弛格子-Boltzm ann模型瞬态辐射常常发生于极短的时间内，在瞬态辐射的模拟中，通常引入无量纲时间来避免过小的时间 步长.将无量纲时间T心代人方程（1)中，得到时间无量纲形式的辐射传递方程[21]di(r,n",t) +L f f. v/(r)/2m,r)= F(r,/r,f*),⑷dt'公式中F{r,n r,r) = i R s{r,{r,r)-L^i(r,f r,r),(5)公式中:心为介质的参考长度.本文提出的时间无量纲形式的辐射传输方程的格子Boltzmann方程如下/(r + c^*,t*+A t')-/(r,t*)=--^(r.t*)] + A t'X),j(6)公式中:/(r，〇为分布函数;M为变换矩阵;S = 士叫U a,…人）为松弛参数矩阵，平衡函数的表达 式为r i(r,n r x)•跑-改2c?辐射强度可以由平衡函数给出，关系如下/(r，/r，〇=-(7)(8)50东北电力大学学报第41卷L B M方法中采用D m Q n格子模型，对于一维和二维问题，本文分别采用D1Q3和D2Q9模型.对于 D1Q3模型,其格子信息为[c〇,c, ,c2] =e；c = [0 1 -l]c,c [2/3,i =0c s=—,(0： = \ll/6,i = 1,2(9)(10)M0 12 一：-1对于D2Q9模型，其格子信息为(11)M C6,C7，C8] y =.0100 1-1-111-11-1l-l- i.c, (12)「4/9，/ =:0ccs = — 〇jt=,1/9,i =]，2,3,4(13).1/36,i=5,6,7,8111111111)-4-1- 1--122224-2- 2-2- 21111010-01—1一 110-20201-1-11•(14) 00 10-111-1-100 - 20211-1-101- 11-1000000 0001-11-h13从格子Boltzm ann方程到辐射传输方程本节基于扩散尺度=7(4幻2下的Maxwell迭代，不引入任何限制和假设，从多松弛格子- Boltzmarm模型严格推导得出辐射传输方程.这种扩散尺度是针对模型中的无量纲时间步长和空间步长 的尺度.首先，令/8U，r))f，w =(叫,叫，…，叫）'时间无量纲形式的 辐射传递方程(6)可以写成矢量形式f(r + ciA t，,r+A t f)-/e9(r,〇] + A t'wF(r,t m),(15)方程（15)左边应用Taylor展开，其中微分算子D' 矩阵/(r + e,A%,«* +y{Ax)2)~^ (A x)*£)s/(r,t*),s = 1〇,=y(E,dx+ E yd y)p(ydt'),P*^s p\q\’A s d i a M e o y e m…，e8,J…，e^)，(16)(17)(18)第1期刘晓川等：求解辐射传输方程的多松弛格子-Bohzmami 模型51公式中和g 均为非负整数.令m = M •/>〃 = M •/%，将Taylor展开形式代入方程（15)并整理得到00工(A x )sDsm =- S [m - me tf] + y (A x )2FMco ,s= 1其中D ^-M D ^-y CE,SX +E ,3y V (y s r yI'*^sp \ q \E t =ME M'E y =M E M '1 .(21)**jJ、’c o定义算子A = X (4幻]5\方程（19)可重新写为s= 1m =m e " -S 'Lm + y (A x )2FS~'Mu , (22)基于扩散尺度下的Maxwell1221迭代，从m° = m 〜开始，方程（19)经过三次迭代得到：m = m" -S ~'[A x D ' + (A x )2D2 + (A x )3D ,]me ，1 + [A x S ^D 1 + (A x )2S ^,Lf2]2ma ,-+7(4.«)2厂5_|财如 + 0((4x )4) ,(23)根据矢量方程(23)的第零项及各算子作用结果，可以得到辐射传递方程a /(r ,/7",f } +L r H" • V /(r ,/T ,<*) = F (r ,/T ,t *) +0((A x )2) ,(24)dt *至此，我们从多松弛格子-Boltzmann模型出发，基于扩散尺度下的Maxwell迭代，严格推导得出了辐射 传输方程，并且可以从方程(24)理论上得出该模型具有二阶的精度.一般而言，对于对流扩散问题，计 算流体力学等问题的L B 模型，其中的松弛系数与宏观方程中的扩散系数，流体黏性系数等有定量关系. 需要指出的是，根据从多松弛格子-Boltzmann模型严格推导得出辐射传输方程可知，本文提出的多松 弛格子-Boltzmann模型中的松弛参数均是自由的，与其他参数无关.对于一维和二维L B 模型，我们取 如下的松弛参数矩阵(19)(20)S = diag( 1 ,ir,l ) ,(25)S = diag(l,l,l,ir,l,5r,1,1,1) ,(26)对流扩散方程的多松弛L B 模型也采用了同样的处理方法，其中一维模型中的松弛参数，二维模型中 的松弛参数h 和s 5与扩散系数有关，而其他的松弛参数均取1.由于松弛矩阵中的松弛参数有无限种组 合方式，因此出于通用性考虑，我们选择了这种处理方法.同时需要指出的是当松弛参数矩阵中的松弛 系数相同时，多松弛模型退化到单松弛模型，即松弛矩阵中的松弛参数均为V2 结果及分析2.1 具有高斯型发射场的一维无限大平板考虑一充满吸收发射性介质的一维无限大平板内的辐射传递问题，平板内具有一高斯型发射场，该 问题由如下方程控制^+/3l ^e -u -b )2/a 2,z,b e [0,1] ,(27)考虑如下边界条件52东北电力大学学报第41卷I (〇,〇 =f 3-]e -b 2/a \ ^>0,(28)该问题存在解析解形式，其表达式如下/(2〇 =/(0，f )e x p ( —,)| 2 - (^ + 6)}X [erf (|+^)-erf (f +a )l ^>0, (29)考虑方向f = 1. 〇，a =〇• 02,6 = 0. 5，采用L B M 来模拟衰减系数为/3 = 1，10和50 时介质内辐射强度的分布，取1〇〇个格子，无量纲时间步长取土‘ =0.000 1，单松弛模型得到的结果和解析解对比，如图1所示,L B M 得到的辐射强度分布和解析解得到的辐射强度分布吻合地很好.接下来，我们进一步研究一维多松弛模型的稳 〇.〇4定性和精度.为了研究稳定性，我们考虑衰减系数为 10 nT1的情况，取100个格子，研究不同松弛参数下|所允许的最大时间步长.数值解和解析解的相对误| 〇.〇2 差定乂为Er = ^-------------(30)丨稳定性标准为数值解和解析解的相对误差小于 10'表1给出了不同松弛参数下所允许的最大时间 步长，不同参数的最大时间步长得到是根据我们定 义的稳定性标准，然后通过数值实验得到的，可以发现多松弛模型允许的最大时间步长可以随松弛参数调整，尤其当松弛参数小于1时，所允许的时间步长 大于单松弛模型，结果表明相比单松弛模型，多松弛模型可以在更大的时间步长内保持稳定，具有更好 的稳定性•多松弛模型的碰撞过程发生在矩空间，与多个速度分布函数相关联，相比单松弛模型发生在 速度空间的碰撞，多松弛模型本身在稳定性方面展现了很大的优势,数值结果证明了多松弛模型在稳定 性上的优势.此外，表2给出了不同格子数下单松弛和多松弛模型的相对误差，可以看出多松弛模型相 比单松弛模型具有更高的精度.图1衰减系数为卢=1,丨〇和501^时LBM 得到的辐射强度分布和解析解对比表1衰减系数/3 = 1〇 n T 1,100个格子下，单松弛（B G K )和多松她（M R T )模型允许的最大时间步长sr =0• 6sr =0. 8 sr =l.O(BGK)sr = 1 • 2sr = l • 4y m a x 22.618.413.28.2 4.1W ax2.26 e-31 • 84 e-31.32 e-30. 82 e-30.41e-3表2衰减系数/3 = 10n不同格子数下,单松弛（B G K )和多松弛（M R T )模型的相对误差格子数sr =0. 65r =0. 85r = l .2sr = 1.4BGK MRT BGKMRT BGK MRT BGK MRT 100 4.24 e-27.72 e-3 1. 14 e -2 3.09 e-3 4. 14 e-3 1.71 e-3 5.79 e-3 3.02 e-3150 1.82 e-2 3.24 e-3 4.84 e-3 1.25 e-3 1.85 e-37.77 e-4 2.54 e-3 1.35 e-32001.01 e-21.79 e-32.68 e~36.83 e—41.04 e-34. 40 e-41.42 e-37.57 e -42.2受高斯型脉冲照射的一维纯散射介质考虑厚度为1 m 的一维半透明平板介质内的瞬态辐射传输问题.介质为各向同性散射，壁面和 介质温度均为〇K，无发射.介质边界为透明边界，环境为真空.平板介质的衰减系数为1 nT1，右侧边界 无照射，左侧边界受到如下法向平行光人射辐射的照射：第1期刘晓川等：求解辐射传输方程的多松弛格子-Boltoimmi 模型53lp(0,t ) = /〇exp [//(〇 ,(31)公式中:/。

流体动力学的格子boltzmann方法及其具体实现格子Boltzmann方法是以Boltzmann方程为基础的，该方程描述了流体中粒子的运动。

格子Boltzmann方法将模拟的流体区域划分为一个个离散的格子，并在每个格子中表示流体的宏观属性，如密度、速度等。

在每个格子中，通过计算碰撞和分布函数来模拟粒子的运动。

具体实现格子Boltzmann方法的步骤如下：1.离散化：首先，将流体区域离散化为一个个格子。

格子的大小可以根据需要进行调整。

2.分布函数：在每个格子中，引入分布函数来描述粒子的密度和速度。

分布函数是一个概率密度函数，表示在给定位置和速度的条件下，粒子在该位置具有该速度的概率。

3.碰撞模拟：在每个格子中，模拟粒子之间的碰撞。

根据碰撞模型，计算粒子之间的相互作用，并更新分布函数。

4.传输：根据速度和分布函数，计算粒子的传输过程。

传输过程描述了粒子从一个格子到另一个格子的流动。

5.边界条件：在模拟流体区域的边界上，需要设置适当的边界条件。

边界条件可以影响流体的流动模式。

6.时间步进：通过迭代计算，不断更新格子中的分布函数。

每个时间步长都对应着碰撞和传输的过程。

格子Boltzmann方法与其他常用的计算流体力学方法相比具有一些优势：1. 高效性：格子Boltzmann方法使用离散化格子的方式来模拟流体运动，计算量相对较小，能够高效地处理大规模流体问题。

2. 并行性：由于格子Boltzmann方法的计算是在各个格子之间进行的，因此可以方便地实现并行计算，利用多核处理器或分布式计算系统，加速计算速度。

3. 多尺度：格子Boltzmann方法可以在不同的尺度上进行模拟，从宏观的流体行为到微观的分子动力学。

4. 可分析性：格子Boltzmann方法建立在Boltzmann方程的基础上，可以通过对方程的分析来推导流体的宏观行为。

总结而言，格子Boltzmann方法是一种基于离散化格子的流体动力学模拟方法，通过计算碰撞和传输过程来模拟流体的运动。

CHEMICAL INDUSTRY AND ENGINEERING PROGRESS 2016年第35卷第6期·1698·化 工 进 展REV 尺度多孔介质格子Boltzmann 方法的数学模型及应用的研究进展张潇丹1，2，雍玉梅2，李文军3，赵元生4，李媛媛2，杨巧文1，杨超2（1中国矿业大学（北京）化学与环境工程学院，北京 100083；2中国科学院过程工程研究所绿色过程与工程重点实验室，北京 100190；3华北科技学院环境工程学院，河北 廊坊 065201；4中国石油化工研究院渣油加氢实验室，北京 102200）摘要：综述了多孔介质表征体元尺度（REV ）格子Boltzmann 模型的研究进展，根据对多孔介质处理方式主要分为部分反弹模型和阻力模型两类，分析归纳了各类模型的优缺点。

由于阻力模型中渗流的广义格子Boltzmann 方程（GLBE ）的作用力是基于GUO 等的作用力模型，可以准确得到宏观方程，不存在离散误差，且模型的平衡分布函数和作用力项中都包含反应介质特性的孔隙率，因而应用最为广泛。

本文还重点介绍了REV 尺度多孔介质LBE 模型在流动、传热、传质、化学反应及相变等过程中的具体应用，认为REV 尺度多孔介质内的三传一反数学模型中需要加入孔隙尺度因素，在更大工程尺度上应该考虑过程参数的各向异性，展望了REV 尺度多孔介质LBE 模型的发展和应用前景。

关键词：多孔介质；表征体元尺度；格子Boltzmann 方法；流动；传热；传质中图分类号：TQ021.9 文献标志码：A 文章编号：1000–6613（2016）06–1698–15 DOI ：10.16085/j.issn.1000-6613.2016.06.010Models and application of lattice Boltzmann method at REV-scalein porous mediaZHANG Xiaodan 1，2，YONG Yumei 2，LI Wenjun 3，ZHAO Yuansheng 4，LI Yuanyuan 2，YANG Qiaowen 1，YANG Chao 2（1School of Chemical & Environmental Engineering ，China University of Mining & Technology (Beijing)，Beijing100083，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 PetroleumProcessing ，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 be第一作者：张潇丹（1989—），女，硕士研究生，主要从事化学工程数值模拟。