计算程序_计算流体力学_对流方程_有限差分法_Lax格式_迎风格式_FTCS格式.docx

% 一维对流方程迎风格式、Lax格式、FTCS格式差分法计算% 潭花林 清华大学航天航空学院% FTCS格式对于一维对流方程不稳定，最好不用clcclear all % 1.参数定义dx=1;x1=-18;x2=18;x=x1:dx:x2;L1=length(x);% dt=0.5*dx; % 收敛dt=2*dx; % 不收敛t1=0;t2=t1+80*dt;t=t1:dt:t2;L2=length(t);alpha=1;lambda=alpha*dt/dx;geshi=1; % 迎风格式% geshi=2; % Lax格式% geshi=3; % FTCS格式 % 2.显式求解zeta=zeros(L1,L2); for kk=1:3 geshi=kk; for ii=1:L1 if x(ii)0 zeta(ii,1)=1; else if x(ii)=0 zeta(ii,1)=1/2; else if x(ii)0 zeta(ii,1)=0; end end endend if geshi=1 for ii=2:L1 for jj=1:(L2-1) zeta(ii,jj+1)=zeta(ii,jj)-lambda*(zeta(ii,jj)-zeta(ii-1,jj); end zeta(1,jj+1)=zeta(2,jj+1); end zeta1=zeta; else if geshi=2 for ii=2:(L1-1) for jj=1:(L2-1) zeta(ii,jj+1)=(zeta(ii+1,jj)+zeta(ii-1,jj)/2-. lambda/2*(zeta(ii+1,jj)-zeta(ii-1,jj); end zeta(1,jj+1)=zeta(2,jj+1); zeta(L1,jj+1)=zeta(L1,jj)-lambda*(zeta(L1,jj)-zeta(L1-1,jj); end zeta2=zeta;else if geshi=3 for ii=2:(L1-1) for jj=1:(L2-1) zeta(ii,jj+1)=zeta(ii,jj)-lambda/2*(zeta(ii+1,jj)-zeta(ii-1,jj); end zeta(1,jj+1)=zeta(2,jj+1); zeta(L1,jj+1)=zeta(L1,jj)-lambda*(zeta(L1,jj)-zeta(L1-1,jj); end zeta3=zeta; end endendend % 3.绘图 % 3.1 t=0figure(1)n=1;plot(x,zeta1(1:L1,n),-k,x,zeta2(1:L1,n),-.k,x,zeta3(1:L1,n),-k)% 作图% axis equal % 是否要求x、y坐标间距相等% grid on % 是否要求画网格xlabel(x/m),ylabel(t/s) % x,y轴表示的变量含义%text(1,2,f(x) % 图中文字标识legend(迎风格式,Lax格式,FTCS格式) % 不同曲线的线型区分title(t=0时刻的计算结果) % 标题 axis(-18,18,-0.2,1.2) % 3.2 t=10figure(2)n=(10-t(1)/dt;plot(x,zeta1(1:L1,n),-k,x,zeta2(1:L1,n),-.k,x,zeta3(1:L1,n),-k)% 作图% axis equal % 是否要求x、y坐标间距相等% grid on % 是否要求画网格xlabel(x/m),ylabel(t/s) % x,y轴表示的变量含义%text(1,2,f(x) % 图中文字标识legend(迎风格式,Lax格式,FTCS格式) % 不同曲线的线型区分title(t=10s时刻的计算结果) % 标题 % 3.3 t=20figure(3)n=(20-t(1)/dt;plot(x,zeta1(1:L1,n),-k,x,zeta2(1:L1,n),-.k,x,zeta3(1:L1,n),-k)% 作图% axis equal % 是否要求x、y坐标间距相等% grid on % 是否要求画网格xlabel(x/m),ylabel(t/s) % x,y轴表示的变量含义%text(1,2,f(x) % 图中文字标识legend(迎风格式,Lax格式,FTCS格式) % 不同曲线的线型区分title(t=20s时刻的计算结果) % 标题 % 3.4 t=40figure(4)n=(40-t(1)/dt;plot(x,zeta1(1:L1,n),-k,x,zeta2(1:L1,n),-.k,x,zeta3(1:L1,n),-k)% 作图% axis equal % 是否要求x、y坐标间距相等% grid on % 是否要求画网格xlabel(x/m),yl