LL用向前差分格式计算初边值问题.doc

偏 微 分 方 程课 程 设 计学号：0683110姓名：陆莉指导老师：翟方曼2010.01一．题目用向前差分格式计算如下热传导方程的初边值问题已知其精确解为二．理论作为模型，考虑一维热传导方程：…………（1.1）其中是正常数，是给定的连续函数现在考虑第二类初边值问题的差分逼近：初始条件：…………（1.2）边值条件：，，………（1.3）假设和在相应区域光滑，并且在满足相容条件，使上述问题有惟一充分光滑的解1.建立差分格式（1）.区域剖分 取空间步长和时间步长，其中都是正整数用两族平行直线和将矩形域分割成矩形网格，网格节点为以表示网格内点集合，即位于开矩形的网点集合；表示所有位于闭矩形的网点集合；是网格界点集合 其次，用表示定义在网点的函数，（2）.微分方程的离散，建立相应差分格式将方程在节点离散化，， …………（1.4）对充分光滑的解，由Taylor展式：…………（1.5）…………（1.6）…………（1.7）（1.5）移项得：…………（1.8）（1.6）（1.7）相加得：…………（1.9）将（1.8）（1.9）代入（1.4）得：…………（1.10）其中，舍去，得到逼近（1.1）的向前格式差分方程：， ……（1.11）其中，，记则由（1.4） 由（1.11） 显然，截断误差（3）.边界条件在本题中，，，，，2.稳定性分析用傅里叶方法对差分格式进行稳定性分析以表示网比，将（1.11）改写成便于计算的形式： （本题中）以代入，得消去，则知增长因子 由，得即 只需 解得 所以向前差分格式的稳定性条件是三． MATLAB程序取，，则，满足稳定性条件另取，，则，亦满足稳定性条件另取，，则，亦满足稳定性条件format longa=2;l=1;T=1;N=10;M=400;h=l/N;to=T/M;r=(a*to)/h^2;for j=1:N+1 x(j)=(j-1)*h; for k=1:M+1 t(k)=(k-1)*to; u(j,k)=exp(x(j)+2*t(k)); endendu %求解精确解for j=1:N+1 x(j)=(j-1)*h; us(j,1)=exp(x(j));endfor k=1:M+1 t(k)=(k-1)*to; us(1,k)=exp(2*t(k)); us(N+1,k)=exp(1+2*t(k));endfor k=2:M+1 for j=2:N us(j,k)=r*us(j-1,k-1)+(1-2*r)*u(j,k-1)+r*us(j+1,k-1); endendus %求解数值解for k=1:M+1 for j=1:N+1 R(j,k)=abs(u(j,k)-us(j,k)); endendR %计算误差Rmax=max(max(R)) %求误差的最大值图---精确解与数值解的比较：x=0:0.1:1;hold onplot(x,u(:,M+1),'b');plot(x,us(:,M+1),'y');title('t=1，h=1/10，τ=1/400时精确解和数值解的比较')text(0.05,21,'蓝：精确解');text(0.05,20,'黄：数值解');hold off图---取不同步长时的误差比较：x=0:1/10:1;y=0:1/20:1;z=0:1/40:1;hold onplot(x,R(:,M+1),'b');hold offM分别取10，20，40四．表格及图表部分结点处的精确解、数值解和误差绝对值（取，）精确解数值解（3，1）（0.2，0）1.2214031.2214030（3，41）（0.2，0.1）1.4918251.4916251.994660e-004（3，81）（0.2，0.2）1.8221191.8218472.719067e-004（3，121）（0.2，0.3）2.2255412.2252053.359077e-004（3，161）（0.2，0.4）2.7182822.7190674.107891e-004（3，201）（0.2，0.5）3.3201173.3196155.018075e-004（3，241）（0.2，0.6）4.0552004.0545886.129183e-004（3，281）（0.2，0.7）4.9530324.9522847.486214e-004（3，321）（0.2，0.8）6.0496476.0487339.143683e-004（3，361）（0.2，0.9）7.3890567.38793911.168120e-004（3，401）（0.2，1）9.0250139.02364913.640773e-004部分结点处的精确解、数值解和误差绝对值（取，）精确解数值解（5，1）（0.2，0）1.2214031.2214030（5，161）（0.2，0.1）1.4918251.4917754.967727e-005（5，321）（0.2，0.2）1.8221191.8220516.790495e-005（5，481）（0.2，0.3）2.2255412.2254578.393560e-005（5，641）（0.2，0.4）2.7182822.71818010.265644e-005（5，801）（0.2，0.5）3.3201173.31999212.540377e-005（5，961）（0.2，0.6）4.0552004.05504715.317112e-005（5，1121）（0.2，0.7）4.9530324.95284518.708399e-005（5，1281）（0.2，0.8）6.0496476.04941922.850495e-005（5，1441）（0.2，0.9）7.3890567.38877727.909658e-005（5，1601）（0.2，1）9.0250139.02467334.088933e-005取不同步长时数值解的最大误差（）23.469161e-0045.911936e-0041.477843e-004图（a）：时精确解和数值解的比较由此图可以看出，精确解和数值解几乎重合。

局部的放大图：由此图可以明显看出数值解和精确解存在的误差图（b）：，取不同步长时精确解与数值解的误差比较蓝色： 绿色 红色由此图可以看出，步长越小，误差越小五． 结论抛物型方程的有限差分法的步骤大致可以归纳如下： 1.对区域进行网格剖分 2.在离散结点建立相应的差分格式 3.处理初边值条件 4.进行稳定性分析由本题可以总结出，抛物型方程的有限差分法所得的数值解能够较好地逼近方程的精确解，且区域剖分得越细，即步长越小，数值解与精确解的误差就越小，数值解越逼近精确解。