1. 引言
Lattice Boltzmann Method (LBM)是在过去的20多年里逐渐发展起来的一种新的不同于传统计算流体力学(CFD)的计算方法 [1] 。相比于采用有限差分、有限体积等传统的计算方法有很大的不同,LBM基于分子动理论,具有比较清楚的物理背景。传统的计算方法是在求解离散方程时直接对宏观的控制方程进行离散化处理,然后使用某种数值方法进行求解,最后得到宏观的各个物理量。而LBM从微观动力学的角度出发,把宏观物理量当作微观量统计平均的结果 [2] 。
LBM方法是由格子气自动机(Lattice Gas Automata, LGA)发展而来 [3] 。在LBM方法中,流体在介观层次上被流散成流体粒子的同时物理区域也被离散为一系列的粒子(微团),时间被离散成一系列的时间步长。其中,描述流体粒子运动的方程成为Boltzmann方程或者相应的离散形式 [4] 。在LBM方法中,速度、密度等宏观物理量不作为变量,基本变量是细观上的粒子分布函数。LBM方法通过粒子分布函数的演变过程,从不同的角度建立起宏观与微观、连续与离散之间的联系,按照简化以后的动力学模型流体粒子被约束在有限的格线上运动,需要对这些相关粒子的相关特征值做平均才可以获得宏观层次的密度、速度等参数 [5] [6] 。
本文针对低雷诺数下的NACA0012翼型进行了数值模拟,通过数值模拟结果与相关文献参考数据的对比,分析并验证LBM方法在数值模拟低雷诺数下翼型绕流时是一种有效、精准的计算方法。
2. 数值方法和几何模型
2.1. 数值方法
LBM法将粒子速度离散为有限速度空间
,N为离散后速度的种类。速度离散后的玻尔兹曼输运方程为
(1)
式中,
为t时刻r点在α方向上的粒子分布函数,
为其离散速度,
为碰撞算子,
为外力在α方向上的分力。
宏观物理量由下式确定:
(2)
(3)
为求解式(1),需对其时间及空间进行离散,得到完全离散的玻尔兹曼方程,并对其进行滤波:
(4)
求解(1)式中必须同时离散速度空间和物理空间,最大的困难在于速度
的处理,由于宏观量是微观量的积分形式,He [7] 用高斯积分,对速度空间离散,得到如下的方程(2D):
(5)
其中:
(6)
(7)
(8)
因为是在2维空间且具有9个离散速度,因此将这种速度离散记为D2Q9。
其实,LBM最初起源于LGA (Lattice gas automata)。经过Higuera和Jimenez提出的一种将碰撞算子线性化的简化模型 [8] 以及Higuera等进一步提出的增加模型的数值稳定性的强化碰撞算子方法 [9] ,LBM去除了统计噪音并去处理碰撞算子指数复杂性的缺点。由于改进后的LBM方法使用的仍然是Fermi-Dirac平衡态分布函数,LGA的其他缺点仍然存在。为了进一步简化碰撞算子,Chen、Qian提出了用SRT (single relaxation time)过程代替碰撞项的格子BGK (LBGK)模型 [10] [11] [12] ,相比之前提出的模型,碰撞算子不会因粒子种类的增加而变得复杂,主要是由于LBGK模型中碰撞算子的本身会发生改变。
通过对速度空间进行离散,宏观量为:
(9)
通过多尺度展开,假定密度变化很小,可由(5)和(9)得到Navier-Stokes方程:
(10)
(11)
这里粘性系数
,压力
。
经过上述处理后,LBM方法完全克服了格子气自动机的一系列缺点。
边界处理是LBM方法进行运用中一项非常重要的内容,它对数值计算的精度、计算的稳定性及其计算的效率有很大的影响。通过LBM方法进行数模拟计算时,由于边界节点上的部分分布函数是未知的,需要根据已知的宏观边界条件确定边界节点的分布函数,因此处理不同边边界条件时需要不同的边界处理格式。目前,将现有的边界处理格式主要分为启发式格式、动力学格式、外推格式以及复杂边界处理格式等四类。相对于其他三种格式,启发式格式不需要较复杂的数学推导的公式求解,一般只要根据边界条件的物理本质通过微观粒子的运动规则直接确定边界节点上的未知分布函数;动力学格式主要是通过边界宏观物理量或者宏观物理量梯度值的定义,直接求解边界节点上未知分布函数的方程组:
,
(12)
以获得边界节点上待定的分布函数。虽然启发式和动力学格式可以处理一些特定的边界条件,但是存在很大的局限性,并且在过程中需要对边界条件的一些物理性质做出假设,很难推广到更一般的边界上。因此,借鉴传统计算流体力学(CFD)方法中的边界处理方法来构造LBM方法边界条件处理格式提出了外推式格式。复杂格式综合了其他三种边界处理格式,主要用于处理曲线边界和运动边界。四类处理格式具有不同的使用范围,其对模拟结果的精度以及稳定影响也差别很大,因此在模拟计算时需要根据具体的物理问题选择合适的边界处理格式。随着LBM方法的不断发展以及解决问题的需要,最近几年多位学者提出了用以处理复杂边界条件的不同格式,如连续反弹格式(continuous bounce-back scheme) [13] 、体积法(volumetric method) [14] [15] 以及浸润边界(immersed boundary) [16] [17] [18] 。
网格划分方面,非均匀网格、多重网格、无网格技术等网格处理技术的出现延伸了LBM方法的应用范围,使得LBM方法逐渐在更多的领域得到实际的应用 [4] 。
2.2. 几何模型
本文选取了NACA0012翼型作为研究对象,它作为一款常用翼型研究数据较多,因此可以方便的进行对比研究。本文算例计算区域大小为60 c × 40 c,计算区域外围为远场边界条件。为了提高数值计算的精度及其可靠性,本文LBM方法所采用的网格生成方法为无网格生成技术。格点在物理量变化较为剧烈的部分分布较密,而变化平缓的部分分布比较稀疏。远场区域解析度为2.56 m,翼型及其尾迹解析度为0.005 m,如图1所示,仅在交界处使用插值的方法,不仅节省计算资源并且保证了计算的精度 [19] 。
3. 数值模拟
本文用LBM方法对雷诺数为
下NACA0012翼型进行数值模拟,将结果与参考文献 [20] [21] 所提供参考值进行对比分析,说明LBM方法在计算低雷诺数下翼型绕流的有效性。为了达到基于弦长等于500的雷诺数,根据表1设定了模拟参数。由于来流为层流,在翼型边界层不应用壁面函数模型。
图2为NACA0012翼型在一流场下不同解析度对于精度的影响的仿真结果。通过不同解析度下LBM法计算得到的阻力系数与CFL3D参考值的对比,结果表明当解析度为0.005 m时,计算结果与文献参考值可以很好的吻合。在翼型及尾迹解析度为0.005 m时,其升力系数、阻力系数如表2所示。
将LBM方法数值模拟流场中五个垂直于翼型截面取无量纲化的X和Y方向速度分量以及压力系数与CFL3D数据进行对比,如图3所示在弦长(
)的五个不同位置。
![](//html.hanspub.org/file/1-2980078x32_hanspub.png)
Figure 1. Grid distribution of NACA0012 airfoil
图1. NACA0012翼型格点分布图
![](//html.hanspub.org/file/1-2980078x33_hanspub.png)
Figure 2. The influence of resolution on drag coefficient
图2. 解析度对阻力系数的影响
图4为LBM法和CFL3D在X与Y方向无量纲化速度分布对比。通过分析我们可以发现:LBM计算所得X方向上的速度分量与CFL3D参考值十分吻合。但LBM法在X方向上的分量与CFL3D参考值在靠近翼型壁面时存在一些偏差,特别是在弦长位置
。这可能是由于计算过程中靠近壁面的粒子层大小的设置不够精准。
压力系数
被定义为
,其中
为静压。图5为LBM法和CFL3D在翼型弦长五个不
同位置压力系数分布对比。通过对比分析可以得知,在计算压力系数时LBM法计算所得结果相对于CFL3D的参考值总体趋势比较吻合,但是在弦长
与
位置靠近壁面处的存在一定的差异。造成壁面处速度的Y方向差异的原因与速度X方向存在差异的原因一样,主要是靠近壁面附近的粒子层精度的设置问题。
表1. 模拟参数
![](Images/Table_Tmp.jpg)
Table 2. NACA0012 lift coefficient and drag coefficient
表2. NACA0012升力系数、阻力系数
综合上述低雷诺数下的NACA0012翼型绕流数值模拟结果的对比分析,可以发现LBM方在对简单对称翼型进行数值模拟具有很高的精准度,模拟结果与文献参考值可以很好地吻合。同时,这也验证了LBM方法在计算低雷诺数下简单对称翼型是一种有效的计算方法。在计算的过程中,对同一流场采用分区域的划分,可以提高区域格点的分布质量,从而在一定程度上能够提高计算的精度。但是,LBM方法在处理翼型边界时采用阶梯近似的方法会对翼型表面物理量的准确性计算产生一定的误差。
4. 结论
Lattice Boltzmann Method (LBM)是最近发展起来的一种新的计算方法,其所具物理图像清晰、边界条件容易处理以并行性能好的特点,使LBM 在一些传统方法研究受到一定限制的领域取得了开拓性的进展。本文采用LBM方法对低雷诺数下的NACA0012翼型进行了模拟,结果表明LBM方法在计算低雷诺数下的简单对称翼型是一种有效的计算方法。未来的工作将以本文为基础,对高雷诺数、复杂边界条件下的翼型进行数值研究,说明LBM方法在更复杂条件下的应用。