偏微分各差分格式matlab编程.docVIP

  • 25
  • 0
  • 约2.32千字
  • 约 5页
  • 2016-12-26 发布于贵州
  • 举报
(1)迎风格式 function u = peHypbYF(a,dt,n,minx,maxx,M) format long; h = (maxx-minx)/(n-1); if a0 for j=1:(n+M) u0(j) = IniU(minx+(j-M-1)*h); end else for j=1:(n+M) u0(j) = IniU(minx+(j-1)*h); end end u1 = u0; for k=1:M if a0 for i=(k+1):n+M u1(i) = -(dt/h)*a*(u0(i)-u0(i-1))+u0(i); end else for i=1:n+M-k u1(i) = -(dt/h)*a*(u0(i+1)-u0(i))+u0(i); end end u0 = u1; end x=0.01:0.01:n/100; u = u1((M+1):(M+n)); plot(x,u,r); axis([0 1 -0.5 1.5]); hold on y=1.*(x=0.5)+0.*(x0.5); plot(x,y,g) title(ó-·???ê?) format short; 运行结果 u=peHypbYF(1,0.005,101,0,1,100) (2)Lax-Fridrichs格式 function u = LaxFridrichs(a,dt,n,minx,maxx,M) format long; h = (maxx-minx)/(n-1); for j=1:(n+2*M) u0(j) = IniU(minx+(j-M-1)*h); end u1 = u0; for k=1:M for i=k+1:n+2*M-k u1(i) = -(dt/h)*a*(u0(i+1)-u0(i-1))/2+(u0(i+1)+u0(i-1))/2; end u0 = u1; end x=0.01:0.01:n/100; u = u1((M+1):(M+n)); plot(x,u,r); axis([0 1 -0.5 1.5]); hold on y=1.*(x=0.5)+0.*(x0.5); plot(x,y,g) title(Lax-Fridrichs??ê?) format short; 运行结果 u = LaxFridrichs(1,0.005,101,0,1,100) (3)Lax-Wendroff格式 function u = LaxWendroff(a,dt,n,minx,maxx,M) format long; h = (maxx-minx)/(n-1); for j=1:(n+2*M) u0(j) = IniU(minx+(j-M-1)*h); end u1 = u0; for k=1:M for i=k+1:n+2*M-k u1(i) = (dt/h)^2*a*a*(u0(i+1)-2*u0(i)+u0(i-1))/2- ... (dt/h)*a*(u0(i+1)-u0(i-1))/2+u0(i); end u0 = u1; end x=0.01:0.01:n/100; u = u1((M+1):(M+n)); plot(x,u,r); axis([0 1 -0.5 1.5]); hold on y=1.*(x=0.5)+0.*(x0.5); plot(x,y,g) title(Lax-Wendroff??ê?) format short 运行结果 u= LaxWendroff(1,0.005,101,0,1,100) (4)Beam-Warming格式 function u =BeamWarming(a,dt,n,minx,maxx,M) format long; h = (maxx-minx)/(n-1); for j=1:(n+2*M) u0(j) = IniU(minx+(j-M-1)*h); end u1 = u0; for k=1:M for i=k+2:n+2*M-k u1(i) = -a*(dt/h)*(u0(i)-u0(i-1))-a*(dt/h)*(1-a*(dt/h

文档评论(0)

1亿VIP精品文档

相关文档