一维热扩散方程的格子Boltzmann方法分析.doc

一维热扩散方程的格子Boltzmann方法分析吴国忠，袁兆成，齐晗兵，李栋，刘杰（东北石油大学 土木建筑工程学院，黑龙江省 大庆市 163318，）摘要：针对一维热扩散方程的数学特点，建立了热扩散方程离散速度模型，构造了其平衡态分布函数，采用Chapman-Enskog 展开和多尺度技术，构建了用于求解一维热扩散方程的D1Q3模型，进行了验证性数值实验。实验结果表明，模型的数值解与文献的解析解吻合良好，其两者的误差随网格细化而大幅度减小，从而说明了本文构建的格子Boltzmann模型可用于求解一维热扩散方程。 关键词：一维热扩散方程；格子Boltzmann方法；Chapman-Enskog 展开；多尺度技术；数值模拟中图分类号：O241.82 文献标识码：ALattice Boltzmann Method Analysis of One Dimension Thermal Diffusion EquationWu Guozhong, Yuan Zhaocheng, Qi Hanbing, Li dong, Liu Jie(Northeast Petroleum University, School of Civil and Architecture Engineering, Daqing, 163318，)Abstract：According to the mathematical characteristics of one dimension thermal diffusion equation, the discrete velocity model and an equilibrium distribution function of thermal diffusion equation were established. D1Q3 model was proposed for one dimension thermal diffusion equation using Chapman-Enskog expansion and multiscale technique. Verification experiments were conducted. The results show that the simulation results are consistent with analytic solutions. The absolute errors between the simulation and analytic solutions become smaller with the mesh refining. And the effectiveness of the lattice Boltzmann model to solve one dimension thermal diffusion equation in this paper is verified.Key words：one dimension thermal diffusion equation; lattice Boltzmann method; Chapman-Enskog expansion; multiscale technique; numerical simulation近年来，格子Boltzmann方法（LBM）在求解偏微分方程领域发展很迅速，特别是求解Navier-Stokes方程获得很大成功。由于Boltzmann方法具基金项目：国家自然科学基金 (No. 51274071)。作者简介：吴国忠（1961-），男，教授，博士生导师，研究方向：格子Boltzmann传热应用； 袁兆成（1990-），男，硕士研究生； 齐晗兵（1975-），男，教授，硕士生导师； 李栋（1979-），男，副教授，硕士生导师； 刘杰（1988-），男，硕士研究生。有物理图像清晰、边界条件容易处理以及并行性能好等优点，被广泛应用于微/纳米尺度流1、多孔介质流2、多相多质流3和晶体生长4等许多领域。扩散过程在物理、化学、生物等许多领域中起着很重要的作用。热扩散方程是描述传热过程的一个重要方程，但在复杂的边界和初始条件下，解析求解是很困难的。许多学者利用格子Boltzmann方法求解扩散问题，并取得了很多成果。刘慕仁等人给出了求解一维有源扩散方程的格子Boltzmann模型，确定了局部平衡函数Chapman-Enskog展开的待定系数5，他们还利用格子Boltzmann方法求解了一维对流扩散方程，确定了方法中的粘滞系数与对流系数的关系6。徐世英等人利用浓度分布的Chapman-Enskog展开及多尺度技术，得出了Tyson反应扩散的反应物和催化剂随时间的浓度空间分布值7。本文在前人的研究基础上，提出了格子Boltzmann方法求解一维热扩散方程的有效数值方法。1 热扩散方程的Boltzmann模型一维热扩散方程表达式： （1）式中，T是温度，要求T 0，u为速度，/cp为热扩散系数，为导热系数，为密度，cp为定压比热容。1.1 平衡态分布函数将一维空间离散成均匀的线性格子，每个节点与相邻的节点相连。将速度ca=cea离散为三个方向，ea可取值-1，0，1，分别对应粒子速度方向向左、不动和向右三种情况。为粒子迁移速率，3-速格子模型如图1所示。图1 3-速格子模型示意图Figure 1. Schematic diagram of 3-bit lattice model用fi(x,t)表示沿ea方向的粒子分布函数，表示在时间t时，在位置x处粒子出现的概率，应有fi0，宏观量T定义为： （2）其中：fieq是平衡态分布函数。根据统计物理，分布函数遵守以下格子Boltzmann方程： （3）其中：是弛豫时间，式（3）也称LBGK（Lattice BGK，简称LBGK）方程，为满足稳定性条件，要求0.5。由式（3）可知，只要得到feq(x,t)就可以得出以后时间的分布函数f(x,t)。由LBGK方程恢复宏观方程的关键是选择适当的平衡态分布函数。由统计力学可知，在局域平衡条件下，分布函数可表示为局域平衡量的函数，考虑到ea的对称性，将平衡态分布函数表示为： （4）式中：a0、a1、a2为待定系数。1.2 格子Boltzmann方程对式（3）的左侧项进行空间和演化时间的Taylor级数展开，取其到二级项可得： （5）应用Chapman-Enskog方法对时间系数、空间系数和分布函数进行多尺度展开，可得到 （6）式中，是个小参量。将式（6）代入（5）式，在一阶小量下有 （7）在二阶小量下有： （8）将（7）式和（8）式分别对下标求和得出 （9） （10）由此得到的式（9）、（10）即为时间t1和t2尺度下的密度平衡方程。1.3 宏观热扩散方程将式（4）代入（2）式，可以得到： （11）为了在t1尺度上恢复方程的对流项，令 （12）由（9）、（12）式得到 （13） （14）由式（7）和式（14）得到 （15）将（15）式代入（10）式，得到 （16）将（14）式乘以，（16）式乘以2，两式相加便可回到宏观尺度下的方程（1）。令，则 （17）由式（11）、（13）和（17）可得出待定系数：，。2 数值实验为验证得到的热扩散方程的Boltzmann模型，采用文献6,8中的算例，进行了验证性模拟实验。2.1 验证算例1 算例1中热扩散方程的边界条件为第一类边界条件且为恒温，具体的控制方程：算例1的解析解在文献6中已给出：算例1中热扩散方程的格子Boltzmann模型求解参数：=1.2，=0.2，=0.5，时间步长=0.01，空间步长=0.1，粒子的迁移速率为c=10，总长度l=40，t分别为280、290、300和310。图2给出了不同时刻t下算例1的模拟值与解析解。从图中可以看出，t为280、290、300和310时，模拟值与解析解的最大绝对误差分别为1.31E-04、8.70E-05、7.20E-05和7.00E-06，从而说明数值解与解析解吻合良好。网格数目对模拟结果精度有一定的影响。本文分析网格数目为80、200和400时，对计算精度的影响情况。 图2 不同时刻下算例1的模拟值与解析解的比较 Figure 2. Comparison of simulation results with analytic solutions at different times in the case 1表1为=1.2，=0.2，=0.5，时间步长=0.01，=300，总长度l=40，空间步长分别为0.5，0.2和0.1时，各个节点在不同网格数目下的解析解与数值解及其绝对误差。从表中可以看出，网格数目分别为80、200和400时的最大绝对误差分别为3.02E-03、3.31E-04和7.00E-05。由此可见，随着网格的细化，模拟结果的绝对误差越来越小。表1各节点在不同网格数目下的解析解与模拟值及其绝对误差Table 1. Simulation results，analytic solutions and absolute error of different points under different grid numbersx解析解80网格200网格400网格模拟值绝对误差模拟值绝对误差模拟值绝对误差50.0010220.0009804.19E-050.0010193.66E-060.0010241.88E-06100.0051340.0049601.74E-040.0051191.50E-050.0051427.92E-06150.0182340.0177255.10E-040.0181904.48E-050.0182562.20E-05200.0536500.0524491.20E-030.0535401.10E-040.0536964.64E-05250.1347350.1324762.26E-030.1345132.22E-040.1348057.00E-05300.2803130.2772963.02E-030.2799823.31E-040.2803655.20E-05350.4123750.4110311.34E-030.4121682.07E-040.4123324.30E-052.2 验证算例2 算例2中的热方程的边界条件是随着时间t变化的，具体方程如下：该问题的解析解在文献8中已给出：图3给出了不同时刻t下算例2的模拟值与解析解。从图中可以看出，t分别为1、4、7和10时的最大绝对误差分别为6.21E-04、4.62E-04、4.98E-04和5.09E-04，从而说明数值解与解析解吻合良好。 图3不同时刻下算例2的模拟值与解析解的比较 Figure 3. Comparison of simulation results with analytic solutions at different times in the case 23 结语本文基于格子Boltzmann方法，采用Chapman-Enskog 展开和多尺度技术，建立了求解一维热扩散方程的D1Q3模型，模型的数值解与文献的解析解吻合良好，其两者的误差随网格细化而大幅度减小，从而说明了本文模型可用于求解一维热扩散方程。参考文献1 Raabe D. Overview of the lattice Boltzmann method for nano- and microscale fluid dynamics in materials science and engineeringJ. Modeling and Simulation in Materials Science Engineering, 2004, 12(6):R13-R46.2 Tang G H, Tao W Q, He Y L. Gas slippage effect on microscale porous flow using the lattice Boltzmann methodJ. Physical Review E, 2005, 72(5):056301.3 Grunau D, Chen S, Eggert K. A lattice Boltzmann model for multiphase fluid flows J. Physics of Fluids, 1993, 5(10):2557-2562.4 Miller W, Succi S, Mansutti D. Lattice Boltzmann model for anisotropic liquid-solid phase transitionJ. Physical review letters, 2001, 86(16):3578-3581.