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

PAGE PAGE 3 学 生 实 验 报 告 实验课程名称 偏微分方程数值解 开课实验室 数统学院 学 院 数 统 年级 2013 专业班 信计02班 学 生 姓 名 学 号 开 课 时 间 2015 至 2016 学年第 2 学期 总 成 绩 教师签名 数学与统计学院制 开课学院、实验室： 数统学院 实验时间 ： 2016年 6月20日 实验项目 名 称 一维波动方程的有限差分法 实验项目类型 验证 演示 综合 设计 其他 指导教师 曾芳 成 绩 是 一．实验目的 通过该实验，要求学生掌握求解一维波动方程的有限差分法，并能通过计算机语言编程实现。 二．实验内容 考虑如下的初值问题： (1) 1．在第三部分写出问题（1）三层显格式。 2．根据你写出的差分格式，编写有限差分法程序。将所写程序放到第四部分。 3．取，分别将时刻的数值解画图显示。 4. 该问题的解析解为，将四个时刻的数值解的误差画图显示，对数值结果进行简单的讨论。 三．实验原理、方法（算法）、步骤 1、三层显格式建立 由于题中，，取，故令网比，，，在内网个点处，利用二阶中心差商得到如下格式： （2） 略去误差项得到： （3） 其中，局部截断误差为。 对于初始条件，建立差分格式为： （4） 对于初始条件，利用中心差商，建立差分格式为： （5） 对于边界条件，建立差分格式为： （6） 将差分格式延拓使为内点，代入（3）得到的式子再与（5）联立消去后整理得到： （7） 综上（3）、（4）、（6）、（7）得到三层显格式如下：（局部截断误差为） （8） 其中。 四．实验环境（所用软件、硬件等）及实验数据文件 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:N u(2,j)=0.5*r^2*u(1,j+1)+(1-r^2)*u(1,j)+0.5*r^2*u(1,j-1); end %定义x=0边界上的数值 for k=1:M+1 u(k,1)=0; end %定义x=1边界上的数值 for k=1:M+1 u(k,N+1)=0; end %迭代计算开始,差分格式 for k=2:M for j=2:N u(k+1,j)=r^2*u(k,j+1)+2*(1-r^2)*u(k,j)+r^2*u(k,j-1)-u(k-1,j); end end u(201,:)=zeros(1,11); %计算k=201行的数值解 u2(201,11)=0; for j=2:N u2(201,j)=r^2*u(200,j+1)+2*(1-r^2)*u(200,j)+r^2*u(200,j-1)-u(199,j); end u=u+u2; u=rot90(u,2);%将矩阵u旋转180度赋值于u %作出图像 [x,t]=meshgrid(0:0.1:1,0:0.01: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)精确解的函数图像); xlabel(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:1; plot(x,wucha(51,:),g*-); hold on grid on plot(x,wucha(101,:),ro-); hold on plot(x,wucha(151,:),ks-); hold on plot(x,wucha(201,:)