热格子Boltzmann法分析及应用

第18卷第5期上濂欢手手报(自然科学版)Vol, 18 No. 52012年10月JOURNAL OF SHANGHAI UNIVERSITY NATURAL SCIENCEoct.2012DOI:10.3969/j.issn.1007-2861.2012.05.010热格子 Boltzmann方法分析及应用陈杰,钱跃竑(上海大学上海市应用数学和力学研究所,上海200072)摘要:格子 boltzmann方法( lattice boltzmann method,LBM)是一种基于气体动理论的介观计算方法,其物理背景清晰、边界处理简单,已成功应用于等温(或无热)流动中简要介绍现有的几种热格子 Boltzmann模型,并运用几种热格子模型求解热 Couette流、方腔自然对流等典型算例,对比不同热格子模型的数值稳定性、准确性、模型的计算效率等将两种热格子模型用于多孔介质内的流动与传热问题中,对比热格子模型在处理复杂结构时的数值特性关键词:格子 Boltzmann方法;热格子 Boltzmann方法;多孔介质中图分类号:0351文献标志码:A文章编号:1007-2861(2012)054048907Analysis and Application of Thermal Lattice Boltzmann MethodCHEN Jie, QIAN Yue-hongShanghai Institute of Applied Mathematics and Mechanics, Shanghai University, Shanghai 200072, China)Abstract: Lattice Boltzmann method LBM)is a mesoscale computational method based on the gaskinetic theory. For solving Fourier- Navier-Stokes equations, the thermal lattice model has attracted muchresearch attention. This paper compares several thermal lattice models in terms of accuracy, stability andcomputational efficiency. The thermal flow in pore-scale porous is also studied using different thermallattice modelsKey words: lattice Boltzmann method( LBM): thermal LBM; porous media格子 Boltzmann方法( lattice Boltzmann method直在不断地探索研究热格子 Boltzmann模型,已形BM)是近20年发展成熟起来的一种数值计算方成了一些经过数值验证具有模拟热流动能力的热法LBM基于气体动理论,通过分布函数的演化获LBM610,并应用于多孔介质流动与传热燃烧及化得宏观信息作为一种简单且能处理复杂流动问题学反应流、湍流等问题.本研究简述了不同热格子的有效数值方法12),LBM具有良好的数值稳定性、 boltzmann模型的基本理论,并通过数值分析对比了天然的并行性、简单的边界处理等优点,自出现之日不同热格子 boltzmann模型的计算结果及数值特性,起就被广泛用于多孔介质流3、多相流反应扩散进而用于多孔介质流动传热问题中系统等诸多领域早期的LBM只应用于等温流动或无热流动)的模拟但是基于这种方法具备处理1等温LBM基本原理复杂问题的能力以及解决传热问题的需要,研究者LBM中问抽数之外,无限维的粒中国煤化工CNMHG收稿日期:201140909基金项目:教育部创新团队资助项目(IRT0844);上海市优秀学术带头人资助项目(12XD1402300)通信作者:钱跃就(1964~),男,教授博士生导师博士研究方向为格子Boltzmann方法理论与应用.Emal:qian@shu.edu.cn490上瀋大報(自然科学版)第18卷子速度空间也都被离散成有限的速度序列.在标准平衡态分布函数是 Maxwell分布的截断形式BM模型中,物理空间被离散成正方形(体)格子,ff =AP+Beipu+流体粒子在格点x上碰撞并按离散速度E=[e,e1,…,en-1]迁移到x+e6,格点f(x,)定义为t时2o<(ciai":8)pu,ug + D.c pu u, (5刻在格点x上速度为e的粒子密度满足如下的格式中,A,Bn,D,为待定参数由满足的守恒条件确子 boltzmann方程:定.平衡态包含了速度的三阶项,离散速度也在f(x+e8,t+δ.)-f(x,t)D2Q9的基础上在主坐标轴上增加了4个速度[f"(x,1)-f(x,)]Qian16采用此模型对一维激波管、二维 Rayleigh式中,为平衡态函数,ω为松弛因子.通过简单地 Benard对流进行了模拟,证明了该模型的有效性向平衡态不断趋近的过程代替真实的复杂碰撞,即MSLBM具有良好的物理基础,宏观方程绝对耦BGK( Bhatnagar-Gross-Krook近似,所以此模型也称合,已成功模拟了一些传热现象但只能模拟狭窄的为LBGK模型.平衡态分布函数的选取是LBM的关温度范围和较小的M数,存在稳定性问题限制了键,DnQm系列中均釆用该模型的广泛应用e;·a(e,·a)2.2熵格子 Boltzman方法(ELBM)熵格子 Boltzmann方法考虑了H定理,通过在式中,,为格子声速W为不同速度粒子的权重本守恒约束下最小化波尔兹曼H函数求解平衡态分布研究在数值模拟中均采用D2Q9模型函数,由此得出的正定的分布函数保证了模型的稳宏观密度和速度分别定义为P=∑,=定性和准确性1. rasianakis等将ELBM拓展到热流动问题的求解中,证实了该方法的有效性,本研∑fe:/p究参照此方法H函数定义为2热格子 Boltzmann模型n=∑/(现有的热格子 Boltzmann模型通常可以分为两大类:第一类是流场温度场耦合统一求解的模型如多平衡态分布函数则是在满足守恒约束条件:∑f=速格子 Boltzmann模型( multi- speed LBM, MSLBM)p,∑f!e.=pu,∑fe2=2pT+pn2的情况下,求H熵格子 boltzmann方法( entropic LBM,ELBM);另一类则是对流场与温度场分别求解如被动标量格子函数最小值得到的具体形式详见文献[10]Boltzmann模型( passive scalar LBM, PSLBM)、双分布Prasianakis等2采用在ELBM中加入高阶量的补函数( double-distribution-function,DF)模型,以及其偿算法较大地提高了基于DQ9标准格子的ELBM可他与传统计算流体动力学( computational fluid模拟的温差和Ma数但是模型实施较为复杂dynamics,CFD)结合的混合方法,如混合热格子2.3双分布函数模型Boltzmann方法( hybrid- thermal LBM, HTLBM双分布函数模型,即存在两个分布函数:密度分2.1多速格子 Boltzmann模型( MSLBM布函数和内能(温度或总能)分布函数,其中密度分多速格子 boltzmann模型是等温LBM模型的直布函数用于模拟速度场,而内能(温度或总能)分布接推广,其密度速度内能等均由速度分布函数的函数则用来模拟温度场.温度、内能或总能分布函数各阶速度矩得到Qian1基于等温LBCK模型,提出均通过不同的方式构造但其演化都独立于密度分了D1Q5,DQ13,D3Q21,D3Q25热力学LBGK模布函数型在这些模型中,除了要满足等温模型的守恒条件231被动标量格子 boltzmann模型(SLBM)外,还应满足能量守恒和平衡态热通量为0的条件被动标量格子 Boltzmann模型基于如下原理:在7=叫(+2"忽略压力做中国煤化工兄下温度可以看作是随CNMHG对流扩散方程-1+)+小,()迪于卖包用周分图模教单为改周第5期陈杰,等:热格子 Boltzmann方法分析及应用题:组分1模拟流体的运动;组分2模拟被动的温度式中,S为广义源项包括压力做的功和粘性热耗散场.平衡态密度函数为速度场与温度场的耦合通过在LBM中添加温n g q= n,W1+(e:u-2}、7)度相关的外力项以及在FVM中添加广义源项S来e·l实现.此外,普朗特数、比热容等热物性以及随温度式中,o表示组分,两组分共享速度,p=∑n,变化的输运系数可以实现相应的调节本研究中FVM与LBM采用同一套网格系统,FVM采用绝对∑稳定且具有与LBM相同精度的二阶迎风格式2.3.2内能双分布函数模型second-order upwind scheme, SUS)内能双分布函数模型最早由He等。提出,其PSLBM,DDF以及 HTLBM这类模型的一个关速度场仍用密度分布函数演化模拟,温度场则由内键之处在于流场与温度场之间的耦合,其模型往往能分布函数模拟该模型的基本思想是通过对连续不满足气体完全状态方程温度场对速度场的影响boltzmann方程进行特殊的离散得到等温LBM,如果只是通过施加一个外力来实现如Cuo等针对进行同样的操作则热LBM可以由离散内能的演化 boussinesq方程组,通过在密度分布函数演化方程中方程得到增加一个外力项以实现温度对流场的影响 Filippova根据内能的定义p=/(-)/4,引人内等“基于HBM研究了小M数下高温燃烧,用温度场修正密度场以满足状态方程.能分布函数g(r,,t)=(-)f2,并引人新的碰撞模型,得到内能分布函数满足的演化方程3计算结果及分析8-8q,(8)为了进一步对比各类模型,本研究采用ELBMatPSLBM,内能DDF模型以及 HTLBM,对热 Couette式中,=(5-)·[0,+(5n]然后对演化方流、封闭方腔自然对流和多孔介质内非等温流动等程离散得到可用于数值计算的离散的分布演化方问题进行了模拟对比程,具体的离散过程详见文献[8]3.1热 Couette流模拟相比于 PSLBM,内能DDF的构造更具有物理基考虑两平板间热 Couette流,上平板以速度U向础,并包含了粘性热耗散和可压缩功相比于右运动,下板静止,且上下平板分别保持恒温T,TMSLBM,DF模型具有更好的数值稳定性,P数不且T>T横截面温度廓线的解析形式为受限制,因此被广泛用于各种近似不可压流体流动PrEc与传热问题2 H(102.4混合热格子 Boltzmann模型( HTLBM)式中,H为平板间距离,Pr=wX为普朗特数,X为热HTLBM是指使用LBM解速度场,使用传统扩散系数E=U/[Cn(T-T)]为埃克特数CFD解温度场,并通过一定的方式相互影响这种方热 Couette流中不考虑流体可压缩性的影响,而法利用了LBM能简单处理复杂流动问题的优势以粘性耗散效应明显,因而分别运用ELBM,内能DDF及传统CFD在传热问题上的成熟技术,可以处理一模型和 HTLBM对该问题进行了模拟,网格数均为些仅仅使用传统CFD较难解决的复杂流动传热问64×64.模拟中Re=UHv=20,计算结果如图1所题最初, Lallemand等将多速多松弛模型和有限示固定P=4,E分别为1,10和20的无量纲温度差分法( finite difference method,FDM)相结合,提出廓线,散点为不同方法的计算值,曲线为解析解公式了混合模型,速度场用多松弛LBM求解,温度场采(10).由图可见,三种模型都成功模拟了粘性耗散用FDM求解效应,且与解析解吻合得很好本研究采用有限容积法( finite volume method,本工作进一步研究了三种模型的计算效率问题FM)与BM相结合的混合方法,即果用如下的图2给出f“中国煤花仝發化曲线,可见FVM求解能量守恒方程:ELBM和H模型CNMHGd(pT),0(pxT)。a(kdr3.2封闭力日热比cpxj封闭方腔尺寸为H(正方形边长),左右壁面分492上降大手军报(自然科学版)第18卷式中B为热膨胀系数物性满足 Boussinesq假设,这O HILBM里通过施加外力G=-B(7-7)g实现温度场对速度场的影响.在方腔自然对流中,可压缩效应以及粘性耗散Ec=10效应可忽略不计.从模型分析可以看出, PSLBM在这种情况下与DDF模型类似,而ELBM边界实施较为复杂因此,本研究分别采用不包含粘性耗散040.60.810效应的 PSLBM和 HTLBM对该问题进行了模拟,模拟中Pr=0.71,Ra数分别为10,103和10图3图1热 Couette流温度廓线和图4分别为HTBM在不同Ra数下流动稳定后得Fig 1 Temperature variation of the thermal Couette now到的流线、等温线,与以往的数值及实验结果一致由图3可见,随着Ra数的增大,方腔中心的近似圆a DDF形的涡逐渐变成椭圆形,进而分裂成两个涡.当Ra=▲ HTLBM10时,两个涡分别向左右壁面移动,在中心出现了第三个涡.由图4可见,随着Ra数的增大,竖直的等温线逐渐变得水平,主导的传热机理由导热变为000对流为了进一步定量考核,本研究计算了努塞尔数№u和平均努塞尔数№um,表1给出了热壁面的Numn、最大Mu数Nm及相应位置的yMm、水平中CPU time/s心线上最大速度tm及相应的位置x、垂直中心线上图2热 Couette流温度残差变化曲线最大速度um以及相应的位置y.HTBM和 PSLBMFig 2 Temperature residuals variation of the therma求解的结果与 Barakos等5的基准解一致Couette flow同样,本研究对 HTLBM和 PSLBM的计算效率别保持恒温T,T,且T>T,上下壁面绝热四壁面进行了对比,图5所示为两种方法模拟自然方腔对速度均为无滑移边界方腔内充满均质空气,考虑向流Ra=10°时,速度残差随CPU时间的变化曲线下的重力描述自然对流的无量纲参数Ra数定义为可以明显看出,两种方法中残差均呈现震荡下降趋势,且 HTLBM收敛快于 PSLBM,HLBM残差收敛RTh-T)HPr到10以下时的耗时为 PSLBM的57%(a)Ra104(b)Ra=103中国煤化工方腔自然对流不同Ra数CNMHGFig 3 Predicted streamlines of natural convection第5期陈杰,等:热格子 boltzmann方法分析及应用493(a)Ra=104(b)Ra=10(c)Ra106图4方腔自然对流不同Ra数的等温线Fig 4 Predicted temperature profiles of natural convection表1数值解与基准解对比Table 1 Comparison of numerical results between thermal models and benchmarksRa数模型№ua(y/Hun(y/HPSLBM2.2473.538(0.1410.194(0.824)0.234(0.121)0.234(0.121)Barakos等163.539(0.140.234(0.119)PSLBM4.5127827(0.0750.128(0.854)0.256(0.065)Ra=107.723(0.085)0.134(0.854)0.260(0.065)4.5107636(0.085)0.132(0.859)0.258(0.066)PSLBM17.454(0.033)0.079(0.852)0.261(0.037)Ra=10°HTLBMl7435(0.040)0.081(0.854)0.263(0.040)Barakos等168.80617442(0.037)0.077(0.859)0.262(0.039)动上有明显的优势及较高的计算率.对于多孔介质内流动与传热的问题,以往使用比较广泛的是BM-5PSLBM和内能DDF模型.本研究将 HTLBM用于多孔介质流动与传热分析中,并与 PSLBM进行了对比本研究分析了分形多孔介质中的自然对流,分形结构采用 Sierpinski地毯,依次对分形等级N=2和3的 Sierpinski情况进行了模拟.无量纲控制参数Pr=0.71,Ra数分别为10,103和10°,固体区域温0100020003000400050006000度保持线性温度分布图6为采用 HTLBM计算NCUP time/s2分形结构内自然对流得到的流线图,图7为相应图5方腔自然对流速度残差变化曲线的等温线由图可见,模拟结果与 PSLBM一致,随rg5 Velocity residuals variation of the natural convection Ra数的逐步增大,传热机理由导热主导变化为对流主导图8为N=3,Ra=10°时的流线图及等温线3.3多孔介质非等温流动模拟由图可见,固体的增多明显地抑制了对流作用多孔介质内部结构十分复杂,其流动传热现象同样对V中国煤化工题上和PSLBM也相当复杂格子 boltzmann方法在模拟孔隙内的流进行了对CNM凵左方法模拟N=2体运动时可以方便地使用反弹格式处理复杂流场分形结构时以是及叫以,此 HTLBM耗时为因此,该方法在孔隙尺度模拟多孔介质内部复杂流FLBM的76%,仍具有优势上藻军手报(自然科学版)第18卷(b)Ra=10(c)Ra=10图6多孔介质方腔自然对流流线(N=2)Fig 6 Predicted streamlines of porous cavity(N=2)(a)Ra=104(b)Ra=105图7多孔介质方腔自然对流等温线(N=2Fig 7 Predicted temperature profiles of porous cavity (N=2)(b)等温线图8多孔介质方腔自然对流流线及等温线(N=3)Fig 8 Predicted streamlines and temperature profiles of porous cavity (N=3)4结论MSLBM, ELH中国煤化工型及 HTLBM),并运用不同CNMHG型算例以及多本研究简要介绍了几种热格子 Boltzmann模型孔介质流动传热问题,得到如下结论第5期陈杰,等:热格子 Boltzmann方法分析及应用495[5]李青,徐旭峰,周美莲三维斑图形成的格子Boltzmann方法模拟[J.上海大学学报:自然科学版,HTLBM2007,13(5):516518[6] QIAN Y H. Simulating thermohydrodynamics with latticeBGK models [J]. Journal of Scientific Computing1993,8(3):231-242[7 SHAN X. Simulation of Rayleigh-Benard convection usinga lattice Boltzmann method [J]. Physical Review E1997,55(3):2780-2788[8 HE X, CHEN S, DOOLEN G D. A novel thermal modelCPU time/sfor the lattice Boltzmann method in incompressible图9多孔方腔自然速度残差变化曲线limit [J]. Joumal of Computational Physics, 1998, 146Fig9 Velocity residuals variation of porous cavity(1):282-300.9 GUO Z, ZHENG C, SHI B, et al. Thermal lattice(1)速度场温度场耦合求解的模型还需要进Boltzmann equation for low Mach number nows步发展才能被广泛应用Decoupling model [J]. Physical Review E, 2007, 75(2)相比于 PSLBM和DDF模型, HTLBM在保(3):036704证计算精度的前提下,具有较高的计算效率[10 PRASIANAKIS N I, CHIKATAMALA S S, KARLIN I V,et(3)数值模拟验证了 HTLBM在处理多孔介质al. Entropic lattice Boltzmann method for simulation of复杂结构时可行、有效,且比 PSLBM的效率高thermal flows [J]. Mathematics and Computers inSimulation,2006,72(2):179-183.[11 ANSUMALI $, KARLIN I V, OTTINGER H C. Minimal參考文献:Europhysics Letters, 2003, 63(6): 798-804[1] QIAN Y H, D'HUMIERES D, LALLEMAND P. Lattice [12] PRASLANAKIS N I, KARLIN I V. Lattice BoltzmannBGK models for Navier-Stokes equation [J]method for simulation of compressible flows on standardEurophysics Letters, 1992, 17(6): 479-484lattices[J]. Physical Review E, 2008, 78(1): 016704[2]QIAN Y H, SUCCI S, ORSZAG SA. Recent advance[13] LALLEMAND P, LUO L S. Theory of the latticelattice Boltzmann computing M]//DIETRICH SBoltzmann method: Acoustic and thermal properties inAnnual reviews of computational physics Il, New Jersey:two and three dimensions [J]. Physical Review EWorld Scientific Publishing Company, 1995: 195-2242003,68(3):0367063] ZHAO C Y, DAI L N, TANGG H, et al. Numerical [14] FILLIPPOVA 0, HANELI D. A novel lattice BGKstudy of natural convection in porous media( metalsapproach for low Mach number combustion [J]. Journalusing lattice Boltzmann method (LBMof Computational Physics, 2000, 158(2): 139-160International Journal of Heat and Fluid Flow, 2010, 31 [15] BARAKOS G, MITSOULIS E, ASSIMACOPOULOS DNatural convection flow in a square cavity revisited[4]严永华,石自媛杨帆液滴撞击液膜喷溅过程的LBMLaminar and turbulent models with wall functions [J]模拟[J].上海大学学报:自然科学版,2008,14(4):International Joumal for Numerical Methods in Fluids,3994041994,18(7):695719中国煤化工CNMHG