格子Boltzmann
格子Boltzmann 方法模拟C/C 复合材料颗粒沉积过程罗思璇()Particle Deposition Process Simulation in C/C Composites by Lattice-Boltzmann MethodLuo Sixuan()Abstract: Lattice Boltzmann method is used here to study the particle deposition process on C/C composites surface. This method considered the boudary condition change during particle deposition. Finally, the deposition pattern is obtained. Keywords: LB Method; flow-particle coupling; C/C composites; deposition摘要:本文使用格子Boltzmann 方法研究了固体火箭发动机中C/C 复合材料表面上颗粒的沉积模态。
该方法考虑了沉积过程中边界形貌的变化对流场的影响,最终得到了颗粒在碳纤维表面的沉积形态。
关键词:LB 方法;流固耦合;C/C 复合材料;沉积0 引言C/C 复合材料是目前新材料领域重点研究和开发的一种新型超高温热结构材料,具有密度小,比强度大、热膨胀系数低、热导率高等特点,是理想的航空航天高温材料[1, 2]。
C/C 复合材料在工作过程中其表面流过的工质为高温燃气。
高温燃气中通常带有燃烧产生的固体颗粒,如选用较高比冲的含铝推进剂时会产生一定量的凝聚相(Al2O3颗粒)。
固体颗粒在C/C 复合材料表面的沉积、冲刷及烧蚀会造成材料内型面的破坏,甚至影响气动性能。
本文使用格子Boltzmann 方法模拟C/C 复合材料中碳纤维上颗粒沉积过程及形态。
1模拟流场的格子Boltzmann 模型格子Boltzmann 方法是近二十年来刚发展起来的,一种以“半晶格分离法”为处理方式的新型热量逐级传递数值方法,最初是在研究电磁场中的流动现象时被提出的,并且该方法可以确定流体域、固体域和温度场在边界处的连续性,十分适合针对复杂几何形状流固耦合传热问题的数值分析。
与传统的经典CFD 方法相比,格子波尔兹曼算法具有很多优点。
因而近年来受到国内外学者的广泛关注,并迅速在气固两相流和传热等研究领域得到应用。
格子Boltzmann 方法将流体抽象为微观的虚拟颗粒,通过这些颗粒在规则的网格点上进行碰撞和迁移来达到模拟流场的目的。
分布函数f i (x ,t )表示t 时刻,x 网格点上,速度为c i 流体颗粒的概率密度,流场的宏观量通过对分布函数进行统计而得到。
本文使用D3Q15模型模拟流场,流体宏观密度ρ和动量ρu 计算如下:10Q i i f ρ-==∑,1Q i i i f ρ-==∑u c(1)本文使用BGK 碰撞算子[3],流场演化方程为:eq (,)(,)[(,)(,)]i i i i i f x t t t f x t f x t f x t τ+⋅∆+∆-=-c(2)其中∆t 为时间步长,τ为无量纲松弛时间,eq i f 为平衡态分布函数,在D2Q9模型中如下计算:2eq2222s s s1[1()]22i i ii f c c c ρα⋅⋅=++-c u c u u(3)其中u 为流体粒子速度,ρ为流体密度;αi 是与模型有关的权系数,D2Q9模型中,α0=4/9,αi =1/9(i=1,3,5,7),αi =1/36(i=2,4,6,8);c s 是当地声速,3s c c ,c =∆x/∆t 。
流体粘度和压力计算公式为:22(21)2s s c P c ντρ⎧=-⎪⎨⎪=⎩(4)2描述颗粒运动的CA 概率模型Dupuis 和Chopard 最早提出了LB-LGA 模型[4, 5],并将其运用到雪的沉积以及水底管道下方沙层的冲刷等定性描述;Przekop 等[6, 7]则利用这种模型定性研究单纤维捕集颗粒过程、纳米和微米尺度混合纤维捕集颗粒过程。
但是这些工作均认为颗粒速度取决于流体微团速度,两者之比为某经验常数,无法定量考虑流体对颗粒的曳力以及其他外力对颗粒运动的影响,无法正确描述颗粒与流体微团之间的轨迹滑移、颗粒与流体之间的相互作用,因此只能定性地描述气固两相流动。
本文使用随机布朗力[8, 9]的方法来描述颗粒的随机布朗扩散运动,因此颗粒的受力方程如下(本文只考虑流体曳力与随机布朗力):p pD B pd dtτ-=+=+u u u F F (5) p p dx u dt=(6)其中u p 为颗粒速度,τp 为颗粒弛豫时间尺度,τp =ρp d p 2/(18μ) ,μ为气体动力粘度,F B 为随机布朗力,ς为平均值0、方差为1的高斯随机数,d p 为颗粒直径,k B 为Boltzmann 常数,T 为温度。
通过对方程(5)进行二次积分计算,可以依次得到颗粒的速度和位移如下1p p f B p ppexp()()(1exp())n n ttu u u F τττ+∆∆=⋅-++⋅⋅--(7) 1p p p f f p B p p p()(1exp())((1exp()))n n n t tx x u u u t t F ττττ+∆∆=+---+⋅∆+∆+--⋅⋅⋅(8)上标n 表示当前时刻,n+1表示下一时刻。
由此可得到颗粒在Δt 内的实际位移Δx (=1p p n nx x +-),把Δx 在某网格方向的投影与该方向网格长度的比值作为颗粒沿该网格线运动到相邻格点的概率:)7,5,3,1(),,0max(=⋅∆=i dxx p ii e (9)其中dx 为网格步长。
最终可以确定颗粒的最终格点位置:1p p 11335577n n x x μμμμ+=++++e e e e(10)其中μi 为一个布尔量,取1的概率为p i 。
以图 1为例,p 1>0, p 3>0, p 5=0, p 7=0(因为e i =-e i +4,p i ×p i +4=0)。
此时通过两个均与分布于[0,1]区间的随机数r 1和r 2,可确定颗粒最终位置,*1123p p*1123p p 1*1123p p 3*1123p p 2,,,,,,,,r p r p r p r p tr p r p t r p r p t⎧>>=⎪⎪<>=+∆⎪⎨><=+∆⎪⎪<<=+∆⎪⎩x x x x e x x e x x e 如果如果如果如果 (11)图 1颗粒运动法则3 颗粒粒径燃烧过程中产生的颗粒粒径分布十分复杂,一般可由三个对数正态分布叠加来表示:23p,p pg,p p 21ppg,pg,ln (/)1()exp 2ln 2ln i i i i iN d d n d d σπσ=⎡⎤=-⎢⎥⎢⎥⎣⎦∑(12)其中,N p,i 为颗粒数目浓度,d pg,i 和σpg,i 分别为颗粒几何平均粒径和标准差,表1所示为各参数具体数值。
颗粒的初始数目浓度分布如图2所示。
表1颗粒尺度分布特征参数参数 N p,1 d pg,1 σpg,1 N p,2 d pg,2 σpg,2 N p,3 d pg,3 σpg,3 值5×10140.081.51×10112.02.01×10910.01.54 计算结果碳纤维的增强形式有单向(1维)、双向(2维)及多向。
单向增强可在一个方向上得到高拉伸强度的C/C ;双向织物提高了抗热应力性能和断裂韧性,容易制成大尺寸形状复杂的部件,有广泛的应用基础。
三向及多向编制具有更好的结构完整性和各向同性。
本文选取双向碳纤维作为模拟对象。
双向增强形式的编织方法有:平纹、缎纹等,本文从中截取一个微元,考察颗粒物在纤维上的沉积模态(见图2)。
图3 双向编织形式中的模拟微元计算区域如图3所示:模拟微元位于流场中央,携带颗粒的气体沿+Z 方向进入计算区域。
为模拟整个碳纤维平面,四周采用周期边界条件。
当颗粒与纤维发生接触时认为发生沉积。
图4展示了颗粒在碳纤维表面的沉积情况。
不同粒径的颗粒由于在流场中的受力不同,沉积的位置也不同。
从图中可见,纤维侧面孔隙沉积颗粒较多,同时承受较大的冲刷。
图3 计算区域图4 颗粒沉积模态5 结论及展望本文使用格子Boltzmann两相流模型对于固体火箭发动机中碳纤维处颗粒沉积现象进行了模拟。
格子Boltzmann方法相比传统的CFD方法,能够将固体边界的形状变化对流场的影响实时反映到计算中。
本文再现了颗粒的沉积过程,在下一步的研究中,可通过实验测量燃烧产生颗粒的真实粒径分布,并通过电镜测量材料表面真实形貌,为更进一步的模拟提供基础。
参考文献:[1] 刘建军,李铁虎,郝志彪. 喉衬热环境与碳/碳复合材料的烧蚀. 宇航材料工艺, 2005, 35 (1): 6.[2] 张智,杨杰,嵇阿琳. 高温处理对轴棒法编织碳/碳复合材料热性能的影响. 材料导报, 2012, 26 (20): 54.[3] Qian, Y. H.,D'Humi猫res, D.,Lallemand, P. Lattice BGK Models for Navier-Stokes Equation. Epl, 1992, 17(6BIS): 479.[4] Dupuis, A.,Chopard, B. Lattice Gas Modeling of Scour Formation under Submarine Pipelines. Journal ofComputational Physics, 2002, 178 (1): 161?174.[5] Masselot, A.,Chopard, B. A lattice Boltzmann model for particle transport and deposition. Epl, 1998, 42 (3): 264.[6] Przekop, R.,Moskal, A.,Grado?, L. Lattice-Boltzmann approach for description of the structure of depositedparticulate matter in fibrous filters. Journal of Aerosol Science, 2003, 34 (2): 133?147.[7] Przekop, R.,Grado?, L. Deposition and Filtration of Nanoparticles in the Composites of Nano- and MicrosizedFibers. Aerosol Science & Technology, 2008, volume 42 (6): 483-493.[8] Maze, B.,Tafreshi, H. V.,Wang, Q., et al. A simulation of unsteady-state filtration via nanofiber media at reducedoperating pressures. Journal of Aerosol Science, 2007, 38 (5): 550?571.[9] Kim, M. M.,Zydney, A. L. Effect of electrostatic, hydrodynamic, and Brownian forces on particle trajectories andsieving in normal flow filtration. J Colloid Interface Sci, 2004, 269 (2): 425?431.。
格子boltzmann方法
格子boltzmann方法
格子Boltzmann方法是一种基于格子模型的统计力学方法,用于计算和模拟多体系统的平衡态和非平衡态性质。
它以物质由大量的微观粒子组成的假设为基础,通过在一个分割成小格子的空间中定义离散的状态,并考虑这些粒子之间的相互作用来描述系统的行为。
在格子Boltzmann方法中,将系统中的宏观性质与微观粒子的状态之间建立联系。
通过定义一个格子上的离散状态,如在每个格子上确定粒子是否存在或具有某些属性,并通过考虑粒子之间的相互作用以及它们在不同的状态之间转移的过程,可以模拟出系统的动力学行为。
这种方法常用于模拟气体动力学、流体力学、固体力学等领域。
格子Boltzmann方法的优点在于它能够处理复杂多体系统,并在很大程度上简化了真实系统的描述。
它可以考虑系统中的不均匀性,如存在的物理场的作用,并可以模拟非平衡态及各种传输过程,如热传导、质量传输等。
格子Boltzmann 方法还可以通过调节格子模型的分辨率以及模型参数的选择来适应不同尺度和
条件下的模拟需求。
然而,格子Boltzmann方法也有一些局限性,如对于高密度和高速度流体的模拟需要更细致的离散化格子,会增加计算复杂度。
此外,由于需要离散化描述系统,格子Boltzmann方法在处理连续和非连续性质之间的界面时可能存在困难。
因此,在具体应用时需要综合考虑这些因素,并结合其他技术和方法进行分析和模拟。
偏微分方程求解的一种新颖方法——格子Boltzmann模型
7 6 邻 节 点
大 学 数 学
第2 7卷
) 撞 , 一 个 节 点 上 从 相 邻 节 点 运 动 来 的 粒 子 发 生 碰 撞 , 据 质 量 、 量 和 能 量 守 恒 规 则 改 碰 在 根 动
其 中 r 松弛 时间 尺度 , 制达 到平衡 的速 度 ( 是 控 可根据 需 要 进行 设 置 ) 由于稳 定 性 的原 因 , 过 实 际测 , 经 算 r必须 大于 1e /.
事 实 上 不 同 的 网 格 剖 分 有 着 不 同 的平 衡 分 布 函数 , B 建 立 模 型 的 核 心 问 题 就 是 根 据 不 同 的 网 格 L M
[ 键 词 ] 格 子 B l ma n方 法 ; 衡 态 分 布 函 数 ; Q 关 ot n z 平 D2 9模 型 ; a i — tk s 程 ; 流一 扩 散 方 程 N ve So e 方 r 对 [ 图 分 类 号 ] O2 1 8 中 4 .2 [ 献标识码]A 文 [ 章 编 号 ] 17 —4 4 2 1 ) 30 7 —8 文 6 21 5 (0 1 0 —0 50
在 低 Mah 马赫 ) 的假 设下 ( l c)其 中粒子平 衡态 分布 函数 为 c( 数 I , U《
~ P
[ + 一 ] +
且
C =c 4 /  ̄,
。 /, 一4 9
1 U 一 3 一(2 一 4 1 9, ∞ 一 6 7 8 1 3 , — / 5 一∞ 一∞ — / 6
第2 7卷 第 3期
格子玻尔兹曼方法
格子玻尔兹曼方法格子玻尔兹曼方法(Lattice Boltzmann Method,简称LBM)是一种基于微观粒子动力学的计算流体力学方法,它通过模拟流体微观粒子在格子空间上的运动来描述流体的宏观行为。
相比传统的有限元方法和有限差分方法,格子玻尔兹曼方法具有较好的并行性能和适应性,特别适用于多孔介质流动、复杂边界条件下的流动以及多相流等问题的模拟。
格子玻尔兹曼方法的基本思想是将流体系统离散化为一个个小的流体微团,这些微团在空间网格上运动,并通过碰撞和迁移过程来模拟流体宏观行为。
在每个时间步长内,微团在空间网格上按照一定的规则进行迁移,并在碰撞过程中遵循玻尔兹曼方程,通过碰撞和迁移过程来模拟流体的宏观行为。
通过在空间网格上迁移和碰撞的过程,可以模拟出流体的宏观运动规律,从而实现对流体流动的模拟和计算。
格子玻尔兹曼方法的优势之一是其较好的并行性能。
由于其基于网格的离散化特性,格子玻尔兹曼方法在并行计算上具有天然的优势,能够有效地利用多核、多节点的计算资源,实现对大规模流体问题的高效模拟。
这使得格子玻尔兹曼方法在计算流体力学领域得到了广泛的应用,特别是在大规模流体模拟和高性能计算方面具有很大的优势。
另外,格子玻尔兹曼方法在处理复杂边界条件和多相流问题上也具有一定的优势。
由于其基于微观粒子动力学的特性,格子玻尔兹曼方法能够比较灵活地处理复杂的边界条件,如固体边界、移动边界等,同时也能够较为方便地模拟多相流体的运动,包括气液两相流、多组分流体等,这使得格子玻尔兹曼方法在工程领域的应用具有广阔的前景。
总的来说,格子玻尔兹曼方法作为一种基于微观粒子动力学的计算流体力学方法,具有较好的并行性能和适应性,特别适用于多孔介质流动、复杂边界条件下的流动以及多相流等问题的模拟。
它在大规模流体模拟和高性能计算方面具有很大的优势,同时也能够比较灵活地处理复杂的边界条件和多相流问题,因此在工程领域具有广泛的应用前景。
格子玻尔兹曼方法的发展将为流体力学领域的研究和工程应用带来新的机遇和挑战。
关于多分布格子boltzmann模型的书
一、概述在统计物理学中,格子Boltzmann模型是一种用于研究粒子在晶格上动力学行为的模型。
在正常的Boltzmann统计力学中,粒子的分布是随机的,而多分布格子Boltzmann模型则引入了多个分布函数,用于描述粒子在不同的晶格上的分布情况。
本文将着重介绍多分布格子Boltzmann模型的相关理论和应用。
二、多分布格子Boltzmann模型的基本概念1. 格子Boltzmann模型的基本原理格子Boltzmann模型最早由硅谷大学的研究者提出,其基本原理是将晶格看作是一个离散的空间,粒子在晶格上的位置和动量也是离散的。
而多分布格子Boltzmann模型则是在每一个晶格上引入一个分布函数,用于描述该格子上粒子的分布情况。
2. 多分布格子Boltzmann模型的表达式多分布格子Boltzmann模型的表达式可以写成如下形式:\[ f_i(\mathbf{r},t) =\sum_{j=1}^{n}\alpha_{ijk}\phi_{ik}(\mathbf{r},t)\]其中,\( f_i(\mathbf{r},t) \)表示晶格i上粒子的分布函数,\( \alpha_{ijk}\)为一个系数,\( \phi_{ik}(\mathbf{r},t) \)为关于晶格i 上粒子的分布函数。
通过引入多个分布函数,我们可以更准确地描述粒子在不同晶格上的动力学行为。
3. 多分布格子Boltzmann模型的演化方程多分布格子Boltzmann模型的演化方程可以写成如下形式:\[ \frac{\partial f_i}{\partial t} + \mathbf{v}_i \cdot \nabla f_i = \frac{1}{\tau_i}(f_{i, eq} - f_i) \]其中,\( f_{i, eq} \)为平衡态分布函数,\( \tau_i \)为弛豫时间。
这个方程描述了不同晶格上粒子的分布函数随时间的演化情况,是多分布格子Boltzmann模型的关键之一。
格子boltzmann方法的原理与应用
格子Boltzmann方法的原理与应用1. 原理介绍格子Boltzmann方法(Lattice Boltzmann Method)是一种基于格子空间的流体模拟方法。
它是通过离散化输运方程,以微分方程的形式描述气体或流体的宏观运动行为,通过在格子点上的分布函数进行更新来模拟流体的动态行为。
格子Boltzmann方法的基本原理可以总结为以下几点:1.分布函数:格子Boltzmann方法中,将流场看作是由离散的分布函数表示的,分布函数描述了在各个速度方向上的分布情况。
通过更新分布函数,模拟流体的宏观行为。
2.离散化模型:为了将连续的流场问题转化为离散的问题,格子Boltzmann方法将流场划分为一个个的格子点,每个格子点上都有一个对应的分布函数。
通过对分布函数进行离散化,实现流场的模拟。
3.背离平衡态:格子Boltzmann方法假设流体运动迅速趋于平衡态,即分布函数以指定的速度在各个方向上收敛到平衡分布。
通过在更新分布函数时引入碰撞过程,模拟流体的运动过程。
4.离散速度模型:分布函数描述了流体在各个速度方向上的分布情况,而格子Boltzmann方法中使用的离散速度模型决定了分布函数的更新方式。
常见的离散速度模型有D2Q9、D3Q15等。
2. 应用领域格子Boltzmann方法作为一种计算流体力学方法,已经在各个领域得到了广泛的应用。
以下是一些常见的应用领域:2.1 流体力学模拟格子Boltzmann方法具有良好的可并行性和模拟精度,适用于复杂流体流动的模拟。
它可以用于模拟包括自由表面流动、多相流动、多物理场耦合等在内的各种复杂流体力学问题。
2.2 细胞生物力学研究格子Boltzmann方法在细胞力学研究中也有广泛应用。
通过模拟流体在细胞表面的流动,可以研究细胞运动、变形和介观流的形成机制。
格子Boltzmann方法在细胞生物力学领域的应用已成为一个重要的研究方向。
2.3 多相流模拟格子Boltzmann方法在多相流动模拟中的应用也非常广泛。
格子boltzmann方法
格子boltzmann方法格子玻尔兹曼方法是一种常用的数值计算方法,它主要用于模拟稀薄气体等流体力学问题。
下面我将从方法原理、模拟过程和应用领域三个方面详细介绍格子玻尔兹曼方法。
首先,格子玻尔兹曼方法基于玻尔兹曼方程和格子Boltzmann方程,通过将连续的物理系统离散化为网格系统进行模拟。
网格系统中的每个格子代表一个微观粒子的状态,而碰撞、传输和外部力的作用通过计算和更新这些格子的状态来实现。
该方法主要包含两个步骤:碰撞和传输。
在碰撞过程中,格子中的粒子通过相互作用和碰撞来改变其速度和方向,从而模拟了分子之间的碰撞过程。
在传输过程中,碰撞后的粒子根据流体的速度场进行移动,从而模拟了背景流场对粒子运动的影响。
其次,在格子玻尔兹曼方法中,模拟的过程可以简化为两个部分:演化和碰撞。
在每个时间步长内,系统首先根据粒子速度和位置的信息计算出相应格点上的分布函数,然后通过碰撞步骤更新这些分布函数以模拟粒子之间的碰撞效应。
通过迭代演化和碰撞步骤,系统的宏观行为可以得到。
格子玻尔兹曼方法中最常用的碰撞操作是BGK碰撞算子,它根据粒子的速度和位置信息计算出新的分布函数,并用该新分布函数代替原来的分布函数。
而在传输过程中,粒子通过碰撞后得到的新速度和方向进行移动。
最后,格子玻尔兹曼方法在流体力学领域具有广泛的应用,特别是在稀薄气体流动、微纳尺度流动和多相流等问题中。
由于其适用于模拟分子尺度和介观尺度流动问题,因此在利用普通的Navier-Stokes方程难以模拟的问题中表现出了良好的效果。
此外,格子玻尔兹曼方法还可以用于模拟流动中的热传导问题、气体分子在多孔介质中的传输问题以及颗粒与流体相互作用等多种复杂流动现象。
近年来,随着计算机性能的不断提高,格子玻尔兹曼方法也得到了快速发展,在模拟大规模真实流体问题方面取得了不错的结果。
总结来说,格子玻尔兹曼方法通过将连续的物理系统离散化为网格系统,模拟粒子碰撞和传输过程,实现了对流体力学问题的数值模拟。
KdV方程的格子Boltzmann模型求解
第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.。
溃坝洪水演进的格子Boltzmann模型
溃坝洪水演进的格子Boltzmann模型摘要:采用基于分子动力学理论的格子boltzmann(lb)方法求解二维浅水方程。
建立了计算二维明渠非恒定流的数学模型。
应用该模型对大坝瞬间全溃或局部溃倒进行了数值模拟,并与理论解进行了比较。
结果说明该模型对模拟溃坝洪水波是有效的。
关键词:格子boltzmann方法二维浅水方程溃坝模拟中图分类号:tv122+.4 文献标识码:a 文章编号:0 引言拦河筑坝可以兴利除害。
同时,大坝溃决也会造成严重的后果。
长期以来,溃坝洪水研究始终是一项重要课题。
随着计算机技术和数值分析方法的发展,数值模拟已经成为研究的主要方法和途径。
二维溃坝波的传播通常由浅水方程来描述。
由于二维浅水方程的一般形式比较复杂,本文中使用无摩擦、无粘性、河底无倾斜的理想情况下的二维浅水方程作为求解模型。
lb方法是近几十年发展起来的一种新的数值计算方法。
它具有算法简单,边界条件易处理,适合并行计算等优点。
本文应用该方法建立了求解二维浅水方程的数学模型,并对计算结果进行了动态演示。
1lb建模1.1 控制方程在水静力学假设条件下,二维浅水波方程可以通过对纳维-斯托克方程沿水深方向进行积分得到。
二维浅水波方程可表示为如下张量形式[1]:式中:和满足爱因斯坦求和约定;表示笛卡尔坐标;表示水深;表示时间;表示方向的垂线平均速度;表示河床高程;表示水流密度;表示水的运动粘度;表示河床的切应力在方向上分量;重力加速度;忽略河床坡率源项和摩阻力源项的理想情况下的二维浅水方程可以表示为:本文将针对上述模型进行建模求解。
1.2 lb方程lb方法基于分子动力学理论,它利用在离散的时空网格上演化的lb方程来模拟宏观流动现象。
带有bgk碰撞项的lattice boltzmann 方程如下[2]:式中为粒子分布函数;为局部平衡态分布函数;为弛豫时间,反映了体系趋于平衡态的快慢。
1.3 chapman-enskog展开和多尺度方程作如下形式的chapman-enskog展开:式中表示。
基于格子Boltzmann方法的单颗粒绕流数值模拟
B ang r o s ok 模 型 和 二 阶精 度 的 曲线 边 界 条 件 处 理 方 法 , 实 现 了对 单 颗 粒 绕 流 问 题 的 定 常 及 非 定 常 h taa- s— o ) Gr Kr
流 动过 程 中 涡结 构 的模 拟 。采 用 动 量 交 换 法 分 别 计 算 了 R y od 数 在 0 1 O en ls . ~2 O范 围 内 2 7个 不 同 R y od 数 时 en ls 的 曳力 系数 ,并 将 计 算 结 果 拟 合得 到 基 于 L M 数 值 模 拟 的 曳 力 曲 线 。计 算 结 果 表 明 , L M 在 气 固 两 相 流 的 模 B B 拟 计 算 中具 有 精 确 、可 靠 的优 点 ,使 用 L M 模 拟 计 算 曳 力 曲线 的方 法 经 济 、 易行 ,并 且 可 以克 服 由传 统 实 验 方 B
非定常单颗粒绕流模拟随着reynolds数的增大颗粒后方的对称旋涡不断增长并出现摆动一般当犚犲100对不稳定的对称旋涡开始渐渐分离形成非对称的周期交替脱落的旋转方向相反的旋涡并以比来流速度小得多的速度向后运动即形成karman给出了在犚犲200时大约一个周期内不同时刻的流线模拟结果可以看出流动分离涡街形成和脱落的过程图中的时间均为量不同reynolds数下定常单颗粒绕流流线图fig4streamlinesofsteadyflowpastsingleparticleatdifferent犚犲5犚犲200不同时刻非定常单颗粒绕流流线图fig5streamlinesofunsteadyflowpastsingleparticleatl w r u i g e pa tc e b s d i l to f fo a o nd a s n l r i l a e
格子Boltzmann方法原理及其应用
格子Boltzmann方法原理及其应用摘要在上世纪八十年代后期提出的格子Boltzamnn方法克服了格子气方法的缺点,其本身也在不断的发展之中.格子Boltzamnn方法在流体运动计算方面展现了非凡的风采,成功地模拟了包括均相不可压缩湍流和多孔介质中的多相流动在内的流体动力学问题.但和成熟的流体动力学计算方法相比,特别在工程实际应用上,该方法还有许多值得研究的地方.本文主要介绍工程实际应用时,具体模型的选择问题.首先从理论上对应用最为广泛的几种基本模型进行了详尽的分析和比较.选择了Poiseuille流动,然后从计算精度、数值稳定性和收敛速度这几个方面进行了细致的比较.从理论和实验两个角度验证了D2G9模型的优越性,为工程实际应用上模型的具体选择提供了一定的参考依据.通过研究二阶精确的格子Boltzamnn模型,提出了非牛顿流体.非牛顿流动性是使用幂法则模型实现的.它可以估算出模型的精确程度,同时不会限制这个模型.二阶精度由剪切变稀和剪切增稠液体的幂法则模型参数范围给出.这些结果与Gabbanelli等人的结果相比,精确度更高,并且得到了更快的计算效率.结果表明了格子Boltzamnn方法适用于非牛顿流体模拟.对于实际流动模拟,本文应用二维9速度模型模拟了四种情况的方柱绕流问题.在第一种情况中,单个方柱位于流场中央,给出了流线图,等涡线图,模拟了卡门涡街现象,并计算了升、阻力系数,Strouhal数等参数;在第二种情况中,计算细长矩板截面柱绕流问题,得到了Strouhal数随着矩形长宽不同的比值下的变化情况;在第三种情况中,两个方柱并列位于流场中央,考察了方柱间距对于流场的影响;在第四种情况中,计算了水平来流为剪切流的方柱绕流问题,比较了速度梯度取不同值下流场的变化情况.所有有关力的求解均采用动量转换法.所得结果,包括流线、等涡线、升/阻力系数曲线等均与已有文献的实验或数值结果基本一致,显示LBM方法及其力的求解方法——动量转换法是有效的,能够精确的模拟各流场.其次,我们还引入一种两相耦合机制对D2G9模型进行了修正,从而使之可以正确处理气固两相流中输运相和颗粒相之间的相互作用.随后,我们模拟了后台阶流动,并和传统CFD方法的模拟结果以及修正其他模型的模拟结果进行了验证,得到了令人满意的结论.从一定程度上验证了两相耦合机制的可行性.通过软件模拟获得了水包油、过渡流型和油包水三种流型的典型模拟图.经分析发现:由软件模拟的流型特点和由探针获得的流型特点具有较好的一致性.在本文最后,我们介绍了以经典算例一方腔流为例,对格子Boltzamnn方法的核心代码进行了优化的方法,主要讲述对时间和空间上的优化,优化的程序使计算效率提高数倍.在并行的框架下,核心演化的代码换为优化后的程序,计算效率有大幅度的提高.关键词:格子方法;格子Boltzamnn 方法;格子气自动机;格子Boltzamnn模型.AbstractIn the latter of 80’s,the Lattice Boltzamnn Method(LBM)was introduced mainlyto cope with major drawbacks of its ancestor,the Lattice Gas Automata(LGA).Eversince,it has undergone a number of refinements and extensions which have taken it tothe point where it can successfully compute a number of non trivial flows,raging fromhomogeneous incompressible turbulence to multiphase flows in porous geometries.Yet,when compared with conventional computational fluids dynamics methods,such as finiteelement,finite difference,it is apparent that there is still a way to go before LBM canachieve full engineering status.In this paper,we mostly focus on the choice of the basic LB models in theengineering application fields.Firstly,we expatiate the basic LB models in theory.Then,we simulate the Poiseuille flow with those basic LB models.And wecompare the simulation results from the computation precision、the numerical stabilityand the convergence rate.Finally,we draw a conclusion that the D2G9 model is the bestchoice in the engineering application fields.Simulation of Flow past square cylinder with LB Method.For the simulation of actual flow,we use D2Q9 investigate fourcases of flow past square cylinders in this paper.For case 1,one singlesquare cylinder is located at the center of the channel,we describe thestreamline contour,vortices contours,simulate the Karman vortex,then compute the lift coefficient,drag coefficient,Strouhal numbersetc.For the case 2,simulate the flow past a cylinder of rectangularcross-section;compute the change of Strouhal numbers varying withthe side ratio.For case 3:two square cylinders arranged side by side inthe center of the channel,the flow features at different spacing ratiosare studied.For case 4:we compute the linear shear flow over a squarecylinder,compare the evolution of flow with different velocitygradient.The results of thesimulation including the streamlines,vorticity contours,lift and drag coefficients etc.are agreed with thoseof available literatures,and show that LB method and itsmomentum-exchange method can achieve accurate results and obtainthe reasonable flow in detail.we employ a two-way coupling mechanisms to modify theD2G9 model.With the modified D2G9 model,we can handle with the interactionsbetween carrier phase and dispersed phase in the model.Then,we simulate abackward-facing step model,and the results are compared qualitatively with the result ofthe traditional CFD method and the other modified LB models.Though the comparison,we can see that the two-way coupling mechanisms can handle with the gas-solid twophases flows successfully.Three kinds of flow pattern,which are oil-in-water flow,transitional flow andwater-in-oil flow,have been got by simulation.According to the result of simulation,theoil-water two-phase flow pattern transition boundary model has been got by.By the analysisof simulation,the characteristic of three kinds of flow pattern of vertical oil which has beengot by analysis of the signals is consistent with results by simulation.We take the classical problem-cavity flow as an example and optimize the kerne codes of the LBM. The optimization include two aspects :time and space .The efficiency of the optimized code increased much more .In the parallel frame,the efficiency also increased if the kernel code is taken the optimized code.Key word:1atrice method;1atrice bohzmann method;lattice gas automata;LBM目录第1章概述 11.1研究格子 Boltzamnn方法的意义 11.2 格子 Boltzamnn方法的发展历程 31.2.1孕育阶段 31.2.2 萌芽到成长阶段 31.3 格子 Boltzamnn方法应用概况及优缺点 51.3.1格子Boltzamnn方法应用概况 51.3.2格子Boltzamnn的优缺点 61.4本论文的研究目的 81.5 相关研究的综述与专注情况 8第2章格子Boltzamnn方法介绍 102.1 Boltzamnn方程的产生 102.2细胞自动机(CA) 112.3格子气自动机(LGA) 122.4格子Boltzamnn方法(LBM) 132.5 格子Boltzamnn的基本结构 162.6本章小结 17第3章格子Boltzamnn方法的基本模型比较 183.1 格子 Boltzamnn 方法基本模型概述 183.2 进行常压力梯度驱动的Poiseuille流动模拟比较几种基本模型 23 3.3本章小结 27第4章格子Boltzamnn方法的算法设计 284.1格子Boltzamnn方法的算法实现 284.2格子Boltzamnn方法的高效算法设计 304.2.1优化算法 304.2.2优化实验 324.3 本章小结 34第5章格子Boltzamnn方法的实际应用 355.1二阶精确格子Boltzamnn非牛顿流体的流动模拟 35 5.1.1理论背景 355.1.2方法和计算结果分析 385.1.3 本节小结 405.2 格子Boltzamnn方法的方柱绕流模拟 405.2.1 单个方柱位于流场中央的绕流问题 405.2.2 细长矩形截面住绕流问题 425.2.3 两个并列方柱的绕流问题 445.2.4来流为剪切流的绕流问题 495.3格子Boltzamnn方法模拟气固两相流 515.3.1对气固两相流的模拟模拟对象简介 515.3.2 计算结果分析 545.3.3本节小结 565.4 格子Boltzamnn方法模拟油水两相流软件设计 565.4.1 LBM油水两相流的关键因素选取 575.4.2 软件的设计 605.4.3 本节小结 635.5 简述格子Boltzamnn方法在其他领域中的应用 645.5.1 颗粒悬浮问题的模拟 645.5.2 热导和对流—扩散问题的模拟 645.5.3 偏微分方程的模拟 655.5.4 多相流和多元流的模拟 65结论及展望 67参考文献 68第1章概述1.1研究格子Boltzamnn方法的意义自从二十世纪四十年代出现了第一台电子计算机以来,人们开始进入了电子信息时代.随着高存储、高速度计算机的出现,人们所能解决的问题也越来越广泛,同时所面临的问题也越来越复杂.在对流动现象的研究中,以往人们大部分依靠的是解析方法,但所解决的问题非常有限.而现实生活中所面临的流动问题往往十分复杂,如航空航天器的亚跨超音速飞行、舰船的航行等等,依靠解析的方法来解决这些复杂的流动现象是不可能的.到现今为止,人们对流体运动的研究主要靠实验方法和数值计算方法.实验方法具有直观、结果基本可靠的特点.但也存在较大的缺点:耗费大、周期长,并且结果受实验条件的影响也较大,尤其是如今的航空航天飞行,速度高、飞行条件复杂,用风洞来模拟困难是相当大的.而流体的运动可以由一组偏微分方程描述.在大多数情况下,这些方程(如N-S方程)都是高度非线性的,采用解析的求解方法是不实际也是不可行的.随着大型计算机的出现,使人们可以借助于计算机用数值计算方法来解决复杂的流动问题.因此,在二十世纪六十年代,用数值方法分析求解流动问题的学科——计算流体力学(CFD)逐渐发展起来.伴随着电子计算机的飞速发展以及各种新颖算法的不断出现,CFD已经形成了一门独立的学科,并且在航空航天、船舶、大型能源装置(如核电站)、新型交通工具、海洋工程、环境保护等众多工程技术部门和领域都得到了广泛的应用.随着计算技术的发展、巨型计算机的出现、计算方法的不断改进,计算流体力学在解决流动的理论和工程实际问题中愈加显示出它的巨大作用.目前,计算流体力学已经成为现代计算科学的最有力的推动力之一.在计算流体力学中,传统的数值模拟方法可以分为两大类:(1)从宏观角度出发,基于连续介质假设,采用数值计算方法,求解全位势方程或Euler方程或N-S方程;(2)从微观角度出发,采用分子动力学的方法,对流动进行数值模拟.其中,格子Boltzamnn方法就是典型的一种.格子Boltzamnn方法(Lattice Boltzamnn Method,LBM)1.1.2格子Boltzamnn法(lattice Boltzamnn method)起源于格子气自动机(Lattice Gas Automata,LGA).LGA方法是元胞自动机(Cellular Automata,CA)在流体力学中的具体应用,是空间、时间和速度空间都离散的一个虚拟微观模型,与以连续微分方程为基础的宏观计算流体力学方法有着本质的不同.LGA的微观特性使得它的边界条件非常容易实现,并且计算也很简单.因此,LGA方法非常适于处理边界复杂的问题.更为重要的是,LGA的计算具有局部性和并行性,非常容易在并行机上实现.LGA的出现不但为并行计算提供了许多新思想,而且对并行计算机制造技术产生了重要的影响.但是,LGA方法也有许多不足之处.例如,由于含有随机因素,LGA的计算结果往往包含很大的统计噪声,LGA的宏观方程也不是标准的流体运动宏观方程.格子Boltzamnn方法是为克服LGA方法的一些内在不足而发展起来的一种新方法.LBM不但克服了LGA的缺点,继承了LGA的主要优点,而且还有许多新的优点,如计算量小、计算效率高、编程简单等.LBM的产生与发展,不仅在计算流体力学领域中产生了深远的影响,它所使用的处理方法和观点对其他许多学科也是富有启发性的.格子Boltzamnn法是一种应用非连续介质思想研究宏观物理现象,并可平行运行,求解流体力学问题的新方法.它是由格子气自动机(lattice gas automata,简称LGA)方法发展而来的.该法把流体及其存在的时间、空间完全离散,把流体看成由许多只有质量没有体积的微小粒子组成,所有这些粒子同步地随着离散的时间步长,根据给定碰撞规则在网格点上相互碰撞,并沿网格线在节点之间运动.碰撞规则遵循质量、动量和能量守恒定律.流体运动的宏观特征是由微观流体格子相互碰撞并在整体上表现出来的统计规律.该法是直接从微观模型出发,经过Boole化处理后进行计算,可认为是N-S差分法逼近的一种无限稳定的格式.被广泛应用于复杂几何边界流体流动、多孔介质流、多相流及反应流等.格子气自动机的基本思想是,把计算区域分成许多均匀的正三角形(或正方形)的网格,而那些只有质量无体积的粒子只能在网格点上存在,并沿着网格线在网格间运动.当某一个粒子从某一网格点到邻近的网格点时,有可能和从其他网格点到达该点的粒子相碰撞.根据Pauli不相容原理,在同一时刻同一点上,沿着每一网格线运动方向最多只有一个粒子,流场中的粒子速度不是0(静止)就是1(设格子边长及时间间隔都为1).以三角形网格为例,每一个网格上在某一时刻,其周围的6个网格上粒子沿着网格线聚集到该点,加上该点可能还有一个静止粒子,这样,可能有7个粒子在该点发生碰撞,然后根据碰撞规则再散射出去,演化为新的运动粒子流向各节点的邻居,形成格子气自动机.1986年MeNamaxa和Zaneltti,提出把格子气自动机中的整数运算变成实数运算,建立了格子Boltzamnn 模型,克服了格子气自动机的数值噪声的缺点.后来陈十一和钱跃宏采用了单一时间松弛方法,满足了各项同性,GalIean不变性,并得到了独立于速度的压力项.使格子Boltzamnn模型保留了格子气自动机的优点,克服了其不足,并在理论分析和数值模拟方面都具有很大灵活性,而且程序编制简单,计算效率较高.从格子Boltzamnn方法诞生至今天已有20年,20年间,其在理论和应用研究等方面都取得了迅速发展,并逐渐成为在相关领域研究的国际热点之一,受到国内外众多学者关注.与之传统模拟方法不同,格子Boltzamnn方法基于分子动理论,具有清晰的物理背景.该方法在宏观上是离散方法,微观上是连续方法,因而被称为介观模拟方法.在许多传统模拟方法难以胜任的领域,入微尺度流动与换热、多孔介质、生物流动、磁流体、晶体生长等,格子Boltzamnn方法都可以进行有效的模拟,因此它被用于多种复杂现象的机理研究,推动了相关学科的发展.可以说,格子Boltzamnn方法不仅仅是一种数值模拟方法,而且是一项重要的科学研究手段.此外,格子Boltzamnn方法还具有天生的并行特性,以及边界条件处理简单、程序易于实施等优点.可以预计,随着计算机技术的进一步发展,以及计算方法的逐渐丰富,格子Boltzamnn方法将会取得更多成果,并为科技发展发挥更重要的作用.1.2 格子Boltzamnn方法的发展历程格子Boltzamnn方法自诞生至今年已取得了长足发展,被誉为现代流体力学的一场变革.1.2.1孕育阶段:对格子Boltzamnn方法发展使得了解,得先从格子自动机说起.格子气自动机使更广泛的元胞自动机在流体学中的应用.元胞自动机是一个时间和空间离散的数学模型.20世纪60年代,Broadwell等人首先提出了离散速度模型,用以研究流体中的激波结构.20世纪70年代,为了研究流体的运输性质,法国的Hardy、Pomeau和Pazzis提出了第一个完全离散模型,该模型命名HPP模型.这是历史上的第一个格子气自动机模型.1986年,法国的Frisch、Pomeau和美国的Hasslacher提出具有足够对称的二维正六变形格子气自动机模型,,命名为FHP模型.由于这些方法在还处在一些缺点:(1)有格子气自动机演化方程推导出来的动量方程不满足Gaililei不变形;(2)流体状态方程不仅仅依赖于密度和温度,还与宏观流速有关;(3)破装蒜子具有指数复杂性,对计算量和存储量也有较大要求.因而,我们将这一段格子气自动机的发展过程称作格子Boltzamnn方法的孕育期.1.2.2 萌芽到成长阶段:自1988年底一篇关于格子Boltzamnn方法的论文出现至今,格子Boltzamnn方法从萌芽逐渐成长壮大,并成为目前一大国际研究热点,受到越来越多学者的关注.1988年,McNamra和Zanetti提出把格子气自动机中的Bool运算变成时数运算,格子点上的粒子数不是用整数0或1来表征,而是用实数f来表示系综平均后的局部粒子分布函数,用Boltzamnn方程代替格子气自动机的演化方程,并将该模型用于流体的数值计算.这是最早的格子Boltzamnn模型,从此开启了格子Boltzamnn方法的历史大门.1989年,Higuera和Jimenez提出了一种简化模型:通过引入平衡分布函数,将碰撞算子线性化.该模型不需要碰撞模型,并忽略各自粒子间的碰撞细节,相比于多粒子碰撞模型,容易构造.同年,Higuera等进一步提出了强化碰撞算子方法,以增加模型的数值稳定性.这两模型统成为矩阵模型.经历了上述两类模型,格子Boltzamnn方法消除了统计噪声,克服了碰撞算子指数复杂性,但是由于依然使用Fermi-Dirac平衡态分布函数,格子气自动机的其他缺点仍然存在.1991年,Chen等提出了单松弛时间法,用同一个时间松弛系数来控制不同例子靠近各自平衡态的快慢,进一步简化了碰撞算子;Qian等人在1992年也提出了类似的方法,称之为格子BGK(LBGK)模型.LBGK模型与矩阵模型类似,但与前面两种模型不同的是,当粒子种类数增加时,碰撞算子本身发生生变化,不会变得复杂.至此,格子Boltzamnn方法完全克服了格子气自动机的一系列缺点,并逐渐成熟,成为国际研究的热点.早期的格子Boltzamnn模型只能用于等温不可压缩流动的模拟.但因为存在可压缩效应,会引起一定的误差.为了消除或强敌有可压缩效应引起的误差,许多学者致力于新的格子Boltzamnn模型的研究,并提出了多种等温不可压模型.而后,一些不可压缩热模型成功实现了对有效范围温度变化的热力学和传热学问题的模型.其中,最成功的要数双分布函数模型.他是在密度分布函数的基础上引入了温度分度函数、或内能分布函数、或总能分布函数,并用密度分布函数演化得到速度场,这类模型具有与等温不可压模型相同的数值稳定性,而且可以从根本上解决压缩功和耗热问题.边界处理方面,经历了20年的发展,格子Boltzamnn方法已逐渐发展出适合不同边界条件、不同模型的边界处理格式.网格划分方面,最初的格子Boltzamnn方法是基于正六边形或正四边形的均匀对称网格.由于均匀网格在计算效率、计算精度等方面的不足,从而促进了非均匀网格、多快以及多重网格、无网格等多技术出现.总的来说,这些网格技术延展了格子Boltzamnn方法的应用范围,使得格子Boltzamnn方法主机去年从理论的神殿走向更可能多的实际应用领域.1.3 格子boltzamnn方法应用概况及优缺点1.3.1格子boltzamnn方法应用概况与传统的宏观数值方法相比,具有介观特性的格子Boltzamnn方法其主要优点是物理图像清晰、便捷容易处理以及并行性能好等.因而自诞生之日起,格子Boltzamnn方法就得到了国内外学术界的广泛关注,并寄希望该方法能再注入为尺度流体、多相流、多孔介质内流动与换热、化学反应流等传统法就延受限的领域取得开拓性进展.事实上,在20年的发展过程中,格子Boltzamnn方法的确也已成一个十分活跃极具发展前景的模拟手段.并迅速在微/纳米尺度流、多孔介质流、多相多质流、非牛顿流体、粒子悬隔i浮流、湍流、化学反应流、燃烧问题、磁流体、晶体生长等许多领域得到应用.下面分别以多孔介质流、多相流和非牛顿流体三个方面为例,做较详细说明.由于格子Boltzamnn方法边界条件易于实施,在模拟具有复杂几何构型的问题具有较大的优势,因而这个方向的发展非常迅速.目前,采用格子Boltzamnn方法对多孔介质流进行模拟主要在空隙尺度和代表单元尺度上进行.在孔隙尺度上,可以直接使用格子Boltzamnn方法描述孔隙内的流体流动,多孔介质则当做固体壁面,流体与介质相互作用使用边界处理格式来描述.在多相流方面,由于真实的流动问题常常是多相的,因而对其开展研究具有重要的现实意义.由于格子Boltzamnn方法的介质特性,它可以方便地描述数流动中不同相之间的相互作用,因而在多相流领域具有较好的应用前景.按照设计方法的不用,现有模拟多相流的格子Boltzamnn模型可分为四大类:着色模型、伪势模型、自由模型和其他模型.格子Boltzamnn方法在非牛顿流体领域的应用刚刚起步,主要研究对象是非牛顿幂律流体.Aharonov等最早提出使用矩阵碰撞该算子来计算幂律流问题,即在每一个时步内,调整碰撞算自来该表局部的动力学黏性系数.Boek用该模型模拟了幂律流体在简化多孔介质中模型的流动,模拟结果与达西定律符合良好.最近,Gabbanelli又对上述模型进行了改进,引入分段幂律方程描述剪切率和表现黏度的关系.以上可看出,到目前为止,格子Boltzamnn方法的研究者主要局限在科学界.尽管如此,随着格子Boltzamnn 方法理论体系逐渐完善,以及计算机技术的进一步发展,格子Boltzamnn方法也会走向更加广泛的工业实际应用中.1.3.2格子Boltzamnn的优缺点流体力学的理论描述通常建立在纳维--斯托克斯方程的基础上,作为流体力学的基石,它已处在了一个多世纪.在通常尺度下,|人们对此方程的物理可靠性即准确性并不抱异议.理论上人们一般通过求纳维--斯托克斯方程及其各种简化形式的途径来处理复杂的流体力学问题,现行的计算流体力学研究也主要是围绕着纳维--斯托克斯方程的计算方法展开的.然而,基于其本质上的非线性以及边界条件处理的困难,除少数简单问题外,解析和数值求解纳维--斯托克斯方程都是极具挑战性的任务.除了求解的困难外,作为一种对流体物理的描述,与描述经典力学运动的牛顿运动方程,或与描述量子力学运动的薛定谔方程等原理方程不同,纳维--斯托克斯方程是从更根本的原理性方程出发,在合理地假定某些物理机制可以忽略后,经过统计平均得到的.本质上纳维--斯托克斯方程当然不可能描述那些被忽略了的物理机制带来的宏观现象,比如流体系统中的相变、非牛顿的本构关系以及在分子运动自由程尺度上的物理现象,在这些领域,纳维--斯托克斯方程明显的显示出了他的局限性.从20世纪80年代末开始,一种对于流体力学的全新的理论表相及有效的计算方法初步形成,这就是现在人们通常所谓的格子Boltzamnn方法.关于格子Boltzamnn方法的早期发展,上文已有较全面的综述,在此仅作简单介绍.从历史角度来讲,格子Boltzamnn方法最初是从所谓的格子气模型演化而来的,而后者是一种抽象简化的分子运动数学模型.格子Boltzamnn方法最初的引入有两个主要原因:一是为了降低模型导致的数值噪音;而是能够克服格子气模型里处在的非物理缺陷.可以证明,格子Boltzamnn系统的宏观表象基本满足纳维--斯托克斯方程.从而,人们可以模拟格子Boltzamnn系统地方法来间接地解纳维--斯托克斯方程.标准格子Boltzamnn方程一般用一下的数学表达式描述:式中——粒子分布函数;——碰撞项.用格子玻尔兹曼模型进行流体的数值模拟有一些明显的优越性.如,它的对流(advection)过程是通过常数值速度实现的.这相应的计算是一项极其简单的操作步骤.当适当的格子网格选定后,该过程通常可以用完全平移的方式实现.用计算数学里的常规有限插值语言来讲,它对应于上风插值.但所不同的是其对应的柯郎数(Courant Number)等于1.相比之下,纳维——斯托克思方程的对流项是一个随时空变化的非线性函数.众所周知,对于它的计算不是一项简单的事,并且,数值稳定性的要求迫使人们在实际问题的计算中只能使用比1小得多的柯朗数.在给定空间分辨度的情况下,小柯朗数意味着小时间步长,从而大大延长了计算时间:同时,小柯朗数也增大了数值扩散误差,迫使人们采用更高精度格式或隐式格式.其后果是,或者算法变得极为复杂,并行效率大大降低;或者计算只限制在处理定常流的情况下.事实上,定常流是对流动情况的极大限制.许多重要的流体力学问题,如分离流,即使我们只关心它的时间平均的结果,也是不能用定常流假设来近似的.在此我们也要提一下格子玻尔兹曼方程的另一个本质特性:所有非线性效应在格子玻尔兹曼方法里都包含在碰撞项中,并且是以纯粹局部信息的方式体现的.这进一步发挥了并行计算的长处.所有这些理由意味着格子玻尔兹曼方法是对非定常流动实行大规模并行模拟计算的一种比较优越的方法.相比之下,以流体力学方程(纳维一斯托克思方程或Burnett类型方程)宏观描述为基础的传统计算方法对许多这类问题存存基本困难.除边界条件之外,利用各种封闭性假设推导出的超越纳维一斯托克思的宏观方程直至现今仍存在对其数学规范性的疑问和争议,多相流的计算也存存同样问题.众所周知,流体系统中存在多相的物理机制是分子问的长程作用力,这种机制早已超出了流体力学方程所能描述的物理现象范围.以流体力学方程为基础的多相流计算方法必须依赖额外的模型来模拟流体力学方程本身所不包含的物理现象.除了实际数值结果显示的问题之外,这种方法本质上隐含着严重的基本物理缺陷,这种缺陷集中表现在对相交界面的准确描述上面,即在十分尖锐的相界面附近,纳维一斯托克思方程之类近平衡态的近似表象是有相当疑问的.这也反映在相界面和兀滑动(no—slip)固体边界条件的互斥性上面,为了修补这一缺憾,人们不得不引入各种滑动经验模型.反之,以细观(mesoscopic)为表象基础的格子玻尔兹曼方法可容忍更大的非平衡态程度及更广义的严格边界条件.另外,压力的状态方程在细观表象中是由粒子的相互作用自然得出的,而不用直接输入和处理.在相变情况下,物体的宏观特性将产生不连续性,而对应的微观和细观力学机制并无改变.格子玻尔兹曼方法在模拟多相流上有着广泛的使用.然而,这种为大多数人所熟悉的格子玻尔兹曼方法的理论框架存在本质上的缺陷.由于它运用逆向切普曼一安斯柯格展开的途径来适定平衡态分布函数中的关键参数,以达到复建宏观物理体系的目的,这就使其。
