格点法matlab代码,Matlab实现格子玻尔兹曼方法.pdf

Matlab实现格子玻尔兹曼方法

Matlab 实现格子玻尔兹曼方法

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% 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 CONSTANTS

lx = 250;

ly = 51;

obst_x = lx/5+1; % position of the cylinder; (exact

obst_y = ly/2+1; % y-symmetry is avoided)

obst_r = ly/10+1; % radius of the cylinder

uMax = 0.02; % maximum velocity of Poiseuille inflow

Re = 100; % Reynolds number

nu = uMax * 2.*obst_r / Re; % kinematic viscosity

omega = 1. / (3*nu+1./2.); % relaxation parameter

maxT = 400000; % total number of iterations

tPlot = 5; % cycles

% D2Q9 LATTICE CONSTANTS

t = [4/9, 1/9,1/9,1/9,1/9, 1/36,1/3 6,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);

% INITIAL CONDITION: (rho=0, u=0) ==> fIn(i) = t(i)

fIn = reshape( t' * ones(1,lx*ly), 9, lx, ly);

% MA