《二维波动方程的有限差分法》由会员分享，可在线阅读，更多相关《二维波动方程的有限差分法（9页珍藏版）》请在金锄头文库上搜索。
1、 学 生 实 验 报 告实验课程名称 偏微分方程数值解 开课实验室 数统学院 学 院 数 统 年级 2013 专业班 信计02班 学 生 姓 名 学 号开 课 时 间 2015至2016学年第2学期总 成 绩教师签名数学与统计学院制开课学院、实验室： 数统学院 实验时间： 2016年 6月20日实验工程名称二维波动方程的有限差分法实验工程类型验证演示综合设计其他指导教师曾芳成 绩是一实验目的通过该实验，要求学生掌握求解二维波动方程的有限差分法，并能通过计算机语言编程实现。二实验容考虑如下的初值问题： (1)1在第三局部写出问题1三层显格式。2根据你写出的差分格式，编写有限差分法程序。将所写程序
2、放到第四局部。3取，分别将时刻的数值解画图显示。4. 该问题的解析解为，将四个时刻的数值解的误差画图显示，对数值结果进展简单的讨论。三实验原理、方法算法、步骤网格划分，故，。在网点，利用二阶中心差商，对1建立差分格式： 2整理得到：3其中，网比，局部截断误差为。考虑边界条件，差分格式为：4考虑初始条件，差分格式为： 5考虑初始条件，利用二阶差商近似：6设时刻的点为点，那么满足差分格式2，代入上式得到：7将6得到的结果代入7中，整理得到：8综上2、4、5、8得到三层显格式的差分格式为： 9其中，局部截断误差为。四实验环境所用软件、硬件等及实验数据文件Matlab%二维波动方程数值计算关键：怎么运
3、用i,j,k三个指标建立循环clc;%可以将代码换成函数m文件h=0.1;tau=0.1*h;%定义步长r=tau/h;%网比x,y,t=meshgrid(0:h:1,0:h:1,0:tau:1.4);%空间网格剖分uu=cos(sqrt(2)*pi*t).*sin(pi*x).*sin(pi*y);%准确解计算%第一层网点计算u=sin(pi*x).*sin(pi*y);%初始条件u1=u(:,:,1);%因为此时得到的u为11x11x141，故只取第一层%第二层网点计算for i=2:10 for j=2:10 u(i,j,2)=0.5*r2*(u(i+1,j,1)+u(i-1,j,1)+
4、u(i,j+1,1)+u(i,j-1,1)+(1-2*r2)*u(i,j,1); u(11,:,2)=0;u(:,11,2)=0; endendu2=u(:,:,2);%第3-141层网点计算for k=2:140 for i=2:10 for j=2:10 u(i,j,k+1)=r2*(u(i+1,j,k)+u(i-1,j,k)+u(i,j+1,k)+u(i,j-1,k)+(2-4*r2)*u(i,j,k)-u(i,j,k-1); u(11,:,k+1)=0;u(:,11,k+1)=0; end endend%结果分析与作图%wucha=abs(u-uu);%求绝对误差矩阵11x11x141
5、wucha1=wucha(:,:,11);%计算t=0.1时刻的绝对误差矩阵11x11wucha2=wucha(:,:,51);%计算t=0.5时刻的绝对误差矩阵11x11wucha3=wucha(:,:,101);%计算t=1.0时刻的绝对误差矩阵11x11wucha4=wucha(:,:,141);%计算t=1.4时刻的绝对误差矩阵11x11x0=0:h:1;y0=0:h:1;%误差分析%作t=0.1时刻的绝对误差图subplot(2,2,1);mesh(x0,y0,wucha1);title(t=0.1时刻的绝对误差);xlabel(x变量);ylabel(y变量);zlabel(绝对误
6、差值);%作t=0.5时刻的绝对误差图subplot(2,2,2);mesh(x0,y0,wucha2);title(t=0.5时刻的绝对误差);xlabel(x变量);ylabel(y变量);zlabel(绝对误差值);%作t=1.0时刻的绝对误差图subplot(2,2,3);mesh(x0,y0,wucha3);title(t=1.0时刻的绝对误差);xlabel(x变量);ylabel(y变量);zlabel(绝对误差值);%作t=1.4时刻的绝对误差图subplot(2,2,4);mesh(x0,y0,wucha4);title(t=1.4时刻的绝对误差);xlabel(x变量);y
7、label(y变量);zlabel(绝对误差值);%四个时刻数值解、准确解%作t=0.1、0.5时刻的数值解与准确解subplot(2,2,1);mesh(x0,y0,u(:,:,11);%作t=0.1时刻的数值解title(t=0.1时刻的数值解);xlabel(x变量);ylabel(y变量);zlabel(u值);subplot(2,2,2);mesh(x0,y0,uu(:,:,11);%作t=0.1时刻的准确解title(t=0.1时刻的准确解);xlabel(x变量);ylabel(y变量);zlabel(u值);%作t=0.5时刻的数值解与准确解subplot(2,2,3);mes
8、h(x0,y0,u(:,:,51);%作t=0.5时刻的数值解title(t=0.5时刻的数值解);xlabel(x变量);ylabel(y变量);zlabel(u值);subplot(2,2,4);mesh(x0,y0,uu(:,:,51);%作t=0.5时刻的准确解title(t=0.5时刻的准确解);xlabel(x变量);ylabel(y变量);zlabel(u值);%分别复制粘贴运行%作t=1.0、1.4时刻的数值解与准确解subplot(2,2,1);mesh(x0,y0,u(:,:,101);%作t=1.0时刻的数值解title(t=1.0时刻的数值解);xlabel(x变量);
9、ylabel(y变量);zlabel(u值);subplot(2,2,2);mesh(x0,y0,uu(:,:,101);%作t=1.0时刻的准确解title(t=1.0时刻的准确解);xlabel(x变量);ylabel(y变量);zlabel(u值);%作t=1.4时刻的数值解与准确解subplot(2,2,3);mesh(x0,y0,u(:,:,141);%作t=1.4时刻的数值解title(t=1.4时刻的数值解);xlabel(x变量);ylabel(y变量);zlabel(u值);subplot(2,2,4);mesh(x0,y0,uu(:,:,141);%作t=1.4时刻的准确解
10、title(t=1.4时刻的准确解);xlabel(x变量);ylabel(y变量);zlabel(u值);五实验结果及实例分析1、时刻的数值解与准确解图 图1 t=0.1、0.5时刻的数值解、准确解图2 t=1.0、1.4时刻的数值解、准确解注：上两图为四个时刻的数值解与准确解，三层显格式达二阶收敛，不难看出，收敛效果很好，符合理论。下列图是四个时刻的绝对误差图像，从图中看出，绝对误差较小，且经过计算得到，收敛阶近似于2，正好符合理论值。2、时刻的绝对误差图图3 四个时刻的绝对误差3、四个时刻t=0.1、0.5、1.0、1.4的绝对误差表t=0.1时刻的绝对误差0.0000 0.0000 0
11、.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0001 0.0001 0.0002 0.0002 0.0002 0.0002 0.0002 0.0001 0.0001 0.0000 0.0000 0.0001 0.0003 0.0004 0.0004 0.0005 0.0004 0.0004 0.0003 0.0001 0.0000 0.0000 0.0002 0.0004 0.0005 0.0006 0.0006 0.0006 0.0005 0.0004 0.0002 0.0000 0.0000
12、0.0002 0.0004 0.0006 0.0007 0.0007 0.0007 0.0006 0.0004 0.0002 0.0000 0.0000 0.0002 0.0005 0.0006 0.0007 0.0008 0.0007 0.0006 0.0005 0.0002 0.0000 0.0000 0.0002 0.0004 0.0006 0.0007 0.0007 0.0007 0.0006 0.0004 0.0002 0.0000 0.0000 0.0002 0.0004 0.0005 0.0006 0.0006 0.0006 0.0005 0.0004 0.0002 0.0000 0.0000 0.0001 0.0003 0.0004 0.0004 0.0005 0.0004 0.0004 0.0003 0.0001 0.0000 0.0000 0.0001 0.0001 0.0002 0.0002 0.0002 0.0002 0.0002 0.0001 0.0001 0.0000 0.0000 0.0000
