一维波动方程的有限差分法.doc

1、学 生 实 验 报 告实验课程名称 偏微分方程数值解 开课实验室 数统学院 学 院 数 统 年级 2013 专业班 信计 02 班 学 生 姓 名 学 号 开 课 时 间 2015 至 2016 学年第 2 学期总 成 绩教师签名数学与统计学院制2开课学院、实验室： 数统学院 实验时间 ： 2016 年 6 月 20 日实验项目类型实验项目名 称 一维波动方程的有限差分法 验证 演示 综合 设计 其他指导教师 曾芳 成 绩 是一实验目的通过该实验，要求学生掌握求解一维波动方程的有限差分法，并能通过计算机语言编程实现。二实验内容考虑如下的初值问题：(1)2,0,1,2,sin,0uxttutut

2、1在第三部分写出问题（1 ）三层显格式。2根据你写出的差分格式，编写有限差分法程序。将所写程序放到第四部分。3取 ，分别将 时刻的数值解画图显示。0.1,.hh0.5,1.,2t4. 该问题的解析解为 ，将四个时刻的数值解的误差画图显示，,cosinuxx对数值结果进行简单的讨论。三实验原理、方法（算法） 、步骤1、三层显格式建立由于题中 ， ，取 ，故令网比 ，0.1,.hh0,1,2xt10,2NM0.1rh， ，在 内网个点处，利用二阶中心差商,2jxj kt jkxt得到如下格式：3（2）11112 22 2kkkkjjj jjjuuoohh略去误差项得到：（3）122211kkkkj

3、j jjjrru 其中 ，局部截断误差为 。1,29,9j oh对于初始条件 ，建立差分格式为：0sinux（4）sin,01,2jjxj对于初始条件 ，利用中心差商，建立差分格式为：,t（5）110,=,2jjjuuj即 ，对于边界条件 ，建立差分格式为：0,2tt（6）0,0kN将差分格式延拓使 为内点，代入（3）得到的式子再与（ 5）联立消去 后整理得到：1ju（7）12202011jj jjurru综上（3） 、 （4） 、 （6） 、 （7）得到三层显格式如下：（局部截断误差为 ）2oh（8）1222110122020110,9,19sinsi,kkkkjj jjjjjjj jjkN

4、rr kuxhrru 其中 。0.rh四实验环境（所用软件、硬件等）及实验数据文件Matlab三层显格式程序如下：%一维波动方程，三层显格式求解法h=0.1;tau=0.1*h;r=tau/h;N=1/h;M=2/tau;x=0:h:1;t=0:tau:2;u=sin(pi*x);%计算 t=0 时刻的 u 值u(1,11)=0;for j=2:Nu(2,j)=0.5*r2*u(1,j+1)+(1-r2)*u(1,j)+0.5*r2*u(1,j-1);end6%定义 x=0 边界上的数值for k=1:M+1u(k,1)=0;end%定义 x=1 边界上的数值for k=1:M+1u(k,N+

5、1)=0;end%迭代计算开始,差分格式for k=2:Mfor j=2:Nu(k+1,j)=r2*u(k,j+1)+2*(1-r2)*u(k,j)+r2*u(k,j-1)-u(k-1,j);endend u(201,:)=zeros(1,11);%计算 k=201 行的数值解u2(201,11)=0;for j=2:Nu2(201,j)=r2*u(200,j+1)+2*(1-r2)*u(200,j)+r2*u(200,j-1)-u(199,j);endu=u+u2;u=rot90(u,2);%将矩阵 u 旋转 180 度赋值于 u%作出图像x,t=meshgrid(0:0.1:1,0:0.0

6、1:2);%划分网格%作出数值解的函数图像subplot(2,2,1);mesh(x,t,u);title(u(x,t)数值解的函数图像);xlabel(x 变量);ylabel(t 变量);zlabel(u 值);%作出精确解的函数图像subplot(2,2,2);u1=cos(pi*t).*sin(pi*x);mesh(x,t,u1);title(u(x,t)精确解的函数图像);7xlabel(x 变量);ylabel(t 变量);zlabel(u 值);%作出 t=0.5，1.0，1.5, 2.0 时刻的绝对误差图像subplot(2,2,3);wucha=abs(u-u1);x=0:h

7、:1;plot(x,wucha(51,:),g*-);hold on grid onplot(x,wucha(101,:),ro-);hold onplot(x,wucha(151,:),ks-);hold onplot(x,wucha(201,:),mp-);title(t=0.5，1.0，1.5, 2.0 时刻的绝对误差函数图像);xlabel(x 变量); ylabel(绝对误差值);legend(t=0.5,t=1.0,t=1.5,t=2.0);%作出 t=0.5，1.0，1.5, 2.0 时刻的数值解函数图像subplot(2,2,4);x=0:h:1;plot(x,u(51,:),

8、g*-);hold ongrid onplot(x,u(101,:),ro-);hold onplot(x,u(151,:),ks-);hold onplot(x,u(201,:),mp-);title(t=0.5，1.0，1.5, 2.0 时刻的数值解函数图像);xlabel(x 变量); ylabel(u 值);legend(t=0.5,t=1.0,t=1.5,t=2.0);%当然也可以作出 u(x,t)绝对误差的函数图像%mesh(x,t,wucha);%title(u(x,t)绝对误差的函数图像);%xlabel(x 变量);%ylabel(t 变量);%zlabel(绝对误差值);3

9、五实验结果及实例分析1、u(x,t)在 t=0.5,1.0,1.5,2.0 时刻的数值解、精确解以及绝对误差表 1 u(x,t)在 t=0.5,1.0,1.5,2.0 时刻的数值解时刻 t t=0.5,1.0,1.5,2.0 时刻的数值解t=0.5 0 -0.0059 -0.0113 -0.0155 -0.0182 -0.0192 -0.0182 -0.0155 -0.0113 -0.0059 0 t=1.0 0 -0.3090 -0.5877 -0.8090 -0.9510 -0.9999 -0.9510 -0.8090 -0.5877 -0.3090 0 t=1.5 0 0.0020 0.

10、0038 0.0052 0.0061 0.0064 0.0061 0.0052 0.0038 0.0020 0 t=2.0 0 0.3090 0.5878 0.8090 0.9511 1.0000 0.9511 0.8090 0.5878 0.3090 0 表 2 u(x,t)在 t=0.5,1.0,1.5,2.0 时刻的精确解时刻 t t=0.5,1.0,1.5,2.0 时刻的精确解t=0.5 0 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0 t=1.0 0 -0.3090 -0.5878 -0.8090 -

11、0.9511 -1.0000 -0.9511 -0.8090 -0.5878 -0.3090 0 t=1.5 0 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0 t=2.0 0 0.3090 0.5878 0.8090 0.9511 1.0000 0.9511 0.8090 0.5878 0.3090 0 表 3 u(x,t)在 t=0.5,1.0,1.5,2.0 时刻的绝对误差时刻 t t=0.5,1.0,1.5,2.0 时刻的绝对误差t=0.5 0 0.0059 0.0113 0.0155 0.0182 0.

12、0192 0.0182 0.0155 0.0113 0.0059 0 t=1.0 0 0.0000 0.0000 0.0001 0.0001 0.0001 0.0001 0.0001 0.0000 0.0000 0 t=1.5 0 0.0020 0.0038 0.0052 0.0061 0.0064 0.0061 0.0052 0.0038 0.0020 0 t=2.0 0 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0 说明：在 t=0.5 时刻的绝对误差最大，t=1.5 时刻次之，t=1 与 t=2 时刻的绝对误差均较小，由于 ，该格式稳定，由数值计算得到的矩阵不难看出，数值解符合理论.1rh解。2、u(x,t)在 t=0.5,1.0,1.5,2.0 时刻的数值解、绝对误差函数图像图 1 数值解、精确解以及绝对误差函数图像3说明：上两图为函数的数值解与精确解，下两图为 t=0.5,1.0,1.5,2.0 时刻的数值解、绝对误差函数图像，符合理论解。教师签名年 月 日