一维抛物线型方程数值解法(1)(附图及matlab程序)(共8页).doc

《一维抛物线型方程数值解法(1)(附图及matlab程序)(共8页).doc》由会员分享，可在线阅读，更多相关《一维抛物线型方程数值解法(1)(附图及matlab程序)(共8页).doc（8页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、精选优质文档-倾情为你奉上 一维抛物线偏微分方程数值解法（1）解一维抛物线型方程(理论书籍可以参看孙志忠：偏微分方程数值解法）Ut-Uxx=0, 0x1,0t0) U(x,0)=ex, 0=x=1,U(0,t)=et,U(1,t)=e(1+t), 0t0)%不用解线性方程组，由下一层(时间层）的值就直接得到上一层的值%m,n为x,t方向的网格数，例如（2-0）/0.01=200;%e为误差，p为精确解u=zeros(n+1,m+1);x=0+(0:m)*h1;t=0+(0:n)*h2;for(i=1:n+1) u(i,1)=exp(t(i); u(i,m+1)=exp(1+t(i);endfo

2、r(i=1:m+1) u(1,i)=exp(x(i);endfor(i=1:n+1) for(j=1:m+1) f(i,j)=0; endendr=h2/(h1*h1); %此处r=a*h2/(h1*h1）；a=1 要求r0)%kmax为最大迭代次数%m,n为x,t方向的网格数，例如（2-0）/0.01=200;%e为误差，p为精确解syms temp;u=zeros(n+1,m+1);x=0+(0:m)*h1;t=0+(0:n)*h2;for(i=1:n+1) u(i,1)=exp(t(i); u(i,m+1)=exp(1+t(i);end for(i=1:m+1) u(1,i)=exp(x

3、(i);endfor(i=1:n+1) for(j=1:m+1) f(i,j)=0; endenda=zeros(n,m-1);r=h2/(h1*h1);%此处r=a*h2/(h1*h1）；a=1 要求rkmax) break; end if(max(max(a) surf(x,t,u) shading interp; xlabel(x);ylabel(t);zlabel(u); title(一维抛物线方程 向前欧拉法 数值解); surf(x,t,p) shading interp;xlabel(x);ylabel(t);zlabel(p); title(一维抛物线方程 向前欧拉法 精确解)

4、同理： plot(x,u) xlabel(x);ylabel(u); title(固定时间 改变x u与x 的关系 数值解)p u e x t=pwxywxq(0.1,0.01,10,100); surf(x,t,u)Warning: Axis limits outside float precision, use ZBuffer or Painters instead. Notrendering Warning: Axis limits outside float precision, use ZBuffer or Painters instead. Notrendering Warning

5、: Axis limits outside float precision, use ZBuffer or Painters instead. Notrendering Warning: Axis limits outside float precision, use ZBuffer or Painters instead. Notrendering surf(x,t,e)Warning: Axis limits outside float precision, use ZBuffer or Painters instead. Notrendering Warning: Axis limits

6、 outside float precision, use ZBuffer or Painters instead. Notrendering Warning: Axis limits outside float precision, use ZBuffer or Painters instead. Notrendering 所以空间步长与时间步长需要满足上面所说的关系 继续减小时间步长p u e x t=pwxywxq(0.1,0.001,10,1000)此为欧拉向前差分法，向后差分法请参看下一篇文章：一维抛物线偏微分方程数值解法（2）（附matlab程序及图片）我近期在做这个，有兴趣可以一起学习百度账号：草随风逝专心-专注-专业