史秀波,闫广武
(1. 桂林理工大学 理学院,广西 桂林 541004;2. 吉林大学 数学学院,长春 130012)
格子Boltzmann方法(LBM)作为一种新的数值方法在计算流体力学、 非线性偏微分方程等领域受到广泛关注[1-2]. 闫广武等[3-4]将该方法应用于波传播问题,为研究其他波动问题提供了可选择的途径. 波动方程中的波速通常是一个常数值,通过引入动量变量ρuj(x,t)将波动方程转换为小扰动Euler方程进行求解[5-6]. 本文用格子Boltzmann方法对热波动方程进行模拟. 在该方程中,波速不再是一个常量,而是一个变量,其表达式为
(1)
其中Cs(x)表示波速,是关于x或y的函数.u(x,t)的下一个时间步表达式为
(2)
本文提出热波的格子Boltzmann模型,通过使用Chapman-Enskog展开和多尺度技术,得到了系列格子Boltzmann偏微分方程、 平衡态分布函数的高阶矩及二阶精度宏观热波方程. 数值实验将模型结果与变分迭代法获得的解析解及经典中心差分格式获得的结果进行比较,结果表明,该方法所得结果与经典方法所得结果相符.
选择一维3-bit网格和二维5-bit网格,分布函数fα(x,t)定义为在某节点x上、t时刻、 具有速度eα(α=0,1,…,b)的粒子出现的概率,其中α=0表示静止粒子. 在一维空间中,b=2,粒子速度为eα={0,c,-c};二维空间中,b=4,粒子速度为
(3)
其中c表示速率. 定义宏观量:
(4)
(5)
格子Boltzmann方程表示为
fα(x+eα,t+1)-fα(x,t)=Ωα+ωα,
(6)
选取Knudsen数ε作为数值模拟的时间步长和Chapman-Enskog展开的小参数[7],在该尺度上,方程(6)可写为
fα(x+εeα,t+ε)-fα(x,t)=Ωα+ωα,
(7)
在方程(7)中,假设
ωα(x,t)=ε2φα(x,t).
(8)
运用Chapman-Enskog展开和多尺度技术,并对方程(7)进行Taylor展开,保留余项到O(ε3)的精度,可得不同时间尺度上的系列格子Boltzmann方程[8]:
结合方程(5),(9),并假设
(12)
可得时间尺度t0上的守恒方程
(13)
式(9)+式(10)×ε并对α求和,同时假设
(14)
得到宏观热波方程为
(15)
结合方程(5),(12),(14)易得平衡态分布函数为
其中D表示空间维数(一维空间中D=1;二维空间中D=2).
为了验证模型效果,分别对一维和二维热波问题进行数值模拟. 一维热波问题使用3-bit模型,二维热波问题使用5-bit模型.
例1一维热波方程
Dirichlet边界条件:
u(0,t)=0,u(1,t)=1+sinht;
(18)b
初始条件:
u(x,0)=x,ut(x,0)=x2.
(18)c
变分迭代法获得的精确解[9]为
u=x+x2sinht.
(18)d
选取参数:格子尺寸m=100,Δx=0.01,c=5.0,τ=1.2,t=1. 图1(A)为t=1时LBM数值解和精确解的比较结果;图1(B)为两种结果的相对误差Er=|(u-u*)/u*|,其中:u表示LBM数值解;u*表示精确解. 由图1(B)可见,相对误差在(1×10-3,9×10-3)内,数值解和精确解吻合较好.
图1 一维热波方程LBM数值解和精确解的比较(A)及相对误差曲线(B)Fig.1 Comparison of LBM solution and the exact solution of one-dimensional thermal wave equation (A) and the curves of their relative error (B)
例2二维热波方程
(19)a
Neumann边界条件:
ux(0,y,t)=0,ux(1,y,t)=2sinht,
uy(x,0,t)=0,uy(x,1,t)=2cosht;
(19)b
初始条件:
u(x,y,0)=y2,ut(x,y,0)=0.
(19)c
选取参数:格子尺寸m×n=100×100,Δx=0.01,Δy=Δx,c=5,τ=1.01,ε=Δt=Δx/c,t=1. 图2(A)为t=1时LBM的模拟结果;图2(B)为t=1时经典中心差分格式的数值解,将其作为精确解;图2(C)为两种结果在x=0.4处的相对误差曲线. 由图2(C)可见,误差区域在(0.00,0.05)内,数值解与精确解吻合较好.
图2 二维热波方程LBM模拟结果(A)、 精确解(B)和两种结果的相对误差曲线(C)Fig.2 LBM result (A),the exact solution (B) and the curves of their relative error (C) for of two-dimensional thermal wave equation
综上,本文提出了一个用于热波方程的格子Boltzmann模型,可得如下结论:
1) 不同时间尺度的系列偏微分方程对构建热波格子Boltzmann模型非常重要,通过使用高阶矩得到了平衡态分布函数的表达式;
[1] CHEN Shi-yi,Doolen G D. Lattice Boltzmann Method for Fluid Flow [J]. Annual Review of Fluid Mechanics,1998,30: 329-364.
[2] DUAN Ya-li,KONG Ling-hua,ZHANG Rui. A Lattice Boltzmann Model for the Generalized Burgers-Huxley Equation [J]. Physica A,2012,391(3): 625-632.
[3] YAN Guang-wu. A Lattice Boltzmann Equation for Waves [J]. Journal of Computational Physics,2000,161(1): 61-69.
[4] YAN Guang-wu,DONG Yin-feng. Application of the Lattice Bhatnagar-Gross-Krook Model to the Simulation of Seismic Pressure Wave [J]. Acta Mechanica Sinica,2005,37(2): 238-243. (闫广武,董银峰. 基于格子Bhatnagar-Gross-Krook模型的地震压力波模拟 [J]. 力学学报,2005,37(2): 238-243.)
[5] ZHANG Jian-ying,YAN Guang-wu,SHI Xiu-bo. Lattice Boltzmann Model for Wave Propagation [J]. Physical Review E,2009,80(2): 026706.
[6] SHI Xiu-bo,YAN Guang-wu,ZHANG Jian-ying. A Multi-energy-level Lattice Boltzmann Model for Two-Dimensional Wave Equation [J]. International Journal for Numerical Methods in Fluids,2010,64(2): 148-162.
[7] Chapman S,Cowling T G. The Mathematical Theory of Non-uniform Gas [M]. Cambridge: Cambridge University Press,1939.
[8] ZHANG Jian-ying,YAN Guang-wu. Lattice Boltzmann Method for One and Two-Dimensional Burgers Equation [J]. Physica A,2008,387(19/20): 4771-4786.
[9] Wazwaz A M,Gorguis A. Exact Solutions for Heat-Like and Wave-Like Equations with Variable Coefficients [J]. Applied Mathematics and Computation,2004,149(1): 15-29.