Matlab实现格子玻尔兹曼方法资料下载.pdf
《Matlab实现格子玻尔兹曼方法资料下载.pdf》由会员分享,可在线阅读,更多相关《Matlab实现格子玻尔兹曼方法资料下载.pdf(3页珍藏版)》请在冰豆网上搜索。
ly=51;
obst_x=lx/5+1;
%positionofthecylinder;
(exactobst_y=ly/2+1;
%y-symmetryisavoided)obst_r=ly/10+1;
%radiusofthecylinderuMax=0.02;
%maximumvelocityofPoiseuilleinflowRe=100;
%Reynoldsnumbernu=uMax*2.*obst_r/Re;
%kinematicviscosityomega=1./(3*nu+1./2.);
%relaxationparametermaxT=400000;
%totalnumberofiterationstPlot=5;
%cycles%D2Q9LATTICECONSTANTSt=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).2fIn(i)=t(i)fIn=reshape(t*ones(1,lx*ly),9,lx,ly);
%MAINLOOP(TIMECYCLES)forcycle=1:
maxT%MACROSCOPICVARIABLESrho=sum(fIn);
ux=reshape(.(cx*reshape(fIn,9,lx*ly),1,lx,ly)./rho;
uy=reshape(.(cy*reshape(fIn,9,lx*ly),1,lx,ly)./rho;
%MACROSCOPIC(DIRICHLET)BOUNDARYCONDITIONS%Inlet:
PoiseuilleprofileL=ly-2;
y=col-1.5;
ux(:
1,col)=4*uMax/(L*L)*(y.*L-y.*y);
uy(:
1,col)=0;
rho(:
1,col)=1./(1-ux(:
1,col).*(.sum(fIn(1,3,5,1,col)+.2*sum(fIn(4,7,8,1,col);
%Outlet:
Zerogradientonrho/uxrho(:
lx,col)=rho(:
lx-1,col);
lx,col)=0;
lx,col)=ux(:
%COLLISIONSTEPfori=1:
9cu=3*(cx(i)*ux+cy(i)*uy);
fEq(i,:
:
)=rho.*t(i).*.(1+cu+1/2*(cu.*cu).-3/2*(ux.2+uy.2);
fOut(i,:
)=fIn(i,:
)-.omega.*(fIn(i,:
)-fEq(i,:
);
end%MICROSCOPICBOUNDARYCONDITIONSfori=1:
9%LeftboundaryfOut(i,1,col)=fEq(i,1,col)+.18*t(i)*cx(i)*cy(i)*(fIn(8,1,col)-.fIn(7,1,col)-fEq(8,1,col)+fEq(7,1,col);
%RightboundaryfOut(i,lx,col)=fEq(i,lx,col)+.18*t(i)*cx(i)*cy(i)*(fIn(6,lx,col)-.fIn(9,lx,col)-fEq(6,lx,col)+fEq(9,lx,col);
%BouncebackregionfOut(i,bbRegion)=fIn(opp(i),bbRegion);
end%STREAMINGSTEPfori=1:
9fIn(i,:
)=.circshift(fOut(i,:
),0,cx(i),cy(i);
end%VISUALIZATIONif(mod(cycle,tPlot)=0)u=reshape(sqrt(ux.2+uy.2),lx,ly);
u(bbRegion)=nan;
imagesc(u);
axisequaloff;
drawnowendend