【电场】基于matlab模拟D2Q9模型LBM圆柱绕流

原创已于 2024-12-23 21:40:21 修改 · 1.1k 阅读

25 ·

CC 4.0 BY-SA版权

文章标签：

#matlab #开发语言

于 2024-03-02 17:34:21 首次发布

物理应用专栏收录该内容

102 篇文章

订阅专栏

本文详细介绍了使用D2Q9模型的格子Boltzmann方法在圆柱绕流模拟中的应用，包括LBM的基本原理、D2Q9模型的实现以及其在数值模拟中的效果。作者还探讨了LBM在这一领域的优势和未来应用前景。

✅作者简介：热爱科研的Matlab仿真开发者，修心和技术同步精进，代码获取、论文复现及科研仿真合作可私信。

🍎个人主页：Matlab科研工作室

🍊个人信条：格物致知。

🔥 内容介绍

本文介绍了使用 D2Q9 模型的格子 Boltzmann 方法 (LBM) 模拟圆柱绕流。LBM 是一种基于粒子动力学的计算流体动力学方法，它可以有效地模拟复杂流动现象。D2Q9 模型是一种常见的 LBM 模型，它具有九个离散速度方向。本文首先介绍了 LBM 的基本原理和 D2Q9 模型的具体实现。然后，本文展示了使用 D2Q9 模型模拟圆柱绕流的数值结果。最后，本文总结了模拟结果并讨论了 LBM 在模拟圆柱绕流中的应用前景。

引言

圆柱绕流是一种常见的流动现象，它广泛存在于工程和自然界中。圆柱绕流的特性对许多工程应用至关重要，例如风力涡轮机、桥梁和建筑物的流体力学设计。传统的计算流体动力学 (CFD) 方法，如有限体积法和有限元法，可以用来模拟圆柱绕流。然而，这些方法通常需要大量的计算资源，并且对于复杂几何形状的流动模拟存在困难。

LBM 是一种基于粒子动力学的 CFD 方法，它可以有效地模拟复杂流动现象。LBM 的基本原理是将流体视为由大量粒子组成的。这些粒子在离散的格子空间中移动，并且相互碰撞。粒子的碰撞遵循特定的规则，这些规则由 LBM 模型定义。D2Q9 模型是一种常见的 LBM 模型，它具有九个离散速度方向。D2Q9 模型可以有效地模拟二维不可压缩流动。

LBM 基本原理

LBM 的基本原理是基于玻尔兹曼方程。玻尔兹曼方程描述了流体中粒子的分布函数随时间和空间的变化。LBM 通过求解玻尔兹曼方程的离散形式来模拟流体流动。

D2Q9 模型的玻尔兹曼方程离散形式为：

f_i(x + e_i, t + 1) - f_i(x, t) = - \frac{1}{\tau} (f_i(x, t) - f_i^{eq}(x, t))

f_i^{eq}(x, t) = \omega_i \rho(x, t) [1 + \frac{e_i \cdot u(x, t)}{c_s^2} + \frac{(e_i \cdot u(x, t))^2}{2c_s^4} - \frac{u(x, t) \cdot u(x, t)}{2c_s^2}]

其中：

圆柱绕流模拟

📣 部分代码

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% cylinder.m: Flow around a cyliner, using LBM        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% This program is free software; you can redistribute it and/or% modify it under the terms of the GNU General Public License% as published by the Free Software Foundation; either version 2% of the License, or (at your option) any later version.% This program is distributed in the hope that it will be useful,% but WITHOUT ANY WARRANTY; without even the implied warranty of% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the% GNU General Public License for more details.% You should have received a copy of the GNU General Public % License along with this program; if not, write to the Free % Software Foundation, Inc., 51 Franklin Street, Fifth Floor,% Boston, MA  02110-1301, USA.%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%clear% GENERAL FLOW CONSTANTSlx         = 250;ly         = 51;obst_x = lx/5+1;   % position of the cylinder; (exactobst_y = ly/2+1;   % y-symmetry is avoided)obst_r = ly/10+1;  % radius of the cylinderuMax  = 0.02;      % maximum velocity of Poiseuille inflowRe     = 100;      % Reynolds numbernu    = uMax * 2.*obst_r / Re;   % kinematic viscosityomega  = 1. / (3*nu+1./2.);      % relaxation parametermaxT   = 400000;   % total number of iterationstPlot  = 5;        % cycles% D2Q9 LATTICE CONSTANTSt  = [4/9, 1/9,1/9,1/9,1/9, 1/36,1/36,1/36,1/36];cx = [  0,   1,  0, -1,  0,    1,  -1,  -1,   1];cy = [  0,   0,  1,  0, -1,    1,   1,  -1,  -1];opp = [ 1,   4,  5,  2,  3,    8,   9,   6,   7];col = [2:(ly-1)];[y,x] = meshgrid(1:ly,1:lx);obst = (x-obst_x).^2 + (y-obst_y).^2 <= obst_r.^2;obst(:,[1,ly]) = 1;bbRegion = find(obst);