二维波动方程的有限差分法

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

2、x,y,t)=0,x,y；：1t10,1.412.3.4.在第三部分写出问题（1）三层显格式。根据你写出的差分格式，编写有限差分法程序。将所写程序放到第四部分。取h=0.14=0.1h,分别将t=0.5,1.0,1.4时刻的数值解画图显示。该问题的解析解为U（x,y,t）=cosJ2：tsinnxsinny,将四个时刻的数值解的误差画图显示，对数值结果进行简单的讨论。.实验原理、方法（算法）、步骤网格划分h=0.1,7=0.1h，故N=1.4htk=kT,k=0,1,111,140。在内网点(K,yj,tk)=10,M=140,Xi=ih,yj=jh,i,j=0,1川，10,利用二阶中心差商，

3、对（1）建立差分格式:k1八kUi,j一2。2k-1kkkkkk'Ui,j=Ui1,j-2ui,j.Ui4,j.Ui,j1-2ui,j,Ui,jh2h2(2)整理得到:k12kkUi,j=rUi1,j'Ui4,jUikj1Uikj42-4r2qkj-UiT(3)其中，i,j=1,2,|119k=1,2川|,139,网比r=2=0.1,局部截断误差为ol+h2）。h考虑边界条件ux,y,t=0,x,y）m.,t.0,1.41,差分格式为：_k_k_k_ku0,0=u0,N=uN,0=uN,N=0,k=0,1JH,140(4)考虑初始条件u（x,y,0）=sinnxsinny,差分

4、格式为:0ui,j=sin二xsinr：yj=sin二ihsin二jh,i,j=0,1,l|l,10(5)2考虑初始条件ut（x,y,0）=0,（x,y产（0,1）,利用二阶差冏近似:11Uj-uij,j,j=0,i,j=0,1,111,102.设k=0时刻的点为内点，则满足差分格式（2）,代入上式得到:u；j=r2（u'j+u：+u：j书+5+）十（24r2从u：(6)将（6）得到的结果u：j=5：代入（7）中，整理得到:1ui,j综上（2）、(4)、k1uui,j12000020三-2r31,jujui,j15,j1-2rui,j(5)、(8)得到三层显格式的差分格式为：2kkkk

5、2kkrui1,juy,jui,j1u,j2-4ru,j-u5(8)0uui,jui1ji,j=1,2,|l|,9,k=1,2川,139u0,0=u0,N=uN,0=uN,N-0,k-0,1,|I,140二sin-xisin.：yjisin二ihsin二jh,i,j=0,1,|l,10=22(50书+u0j书+u：j)+(12r2以/，j=0,1,|,10(9)其中=三=0.1,局部截断误差为o（T2+h2）。h四.实验环境（所用软件、硬件等）及实验数据文件Matlab%二维波动方程数值计算（关键：怎么运用i,j,k三个指标建立循环）clc;%可以将代码换成函数m文件h=0.1;tau=0.1

6、*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,故只取第一层%第二层网点计算fori=2:10forj=2:10u(i,j,2)=0.5*rA2*(u(i+1,j,1)+u(i-1,j,1)+u(i,j+1,1)+u(i,j-1,1)+(1-2*rA2)*u(i,j,1)；u(11,:

7、,2)=0;u(:,11,2)=0;endendu2=u(:,:,2);%S3-141层网点计算fork=2:140fori=2:10forj=2:10u(i,j,k+1)=rA2*(u(i+1,j,k)+u(i-1,j,k)+u(i,j+1,k)+u(i,j-1,k)+(2-4*rA2)*u(i,j,k)-u(i,j,k-1);u(11,:,k+1)=0;u(:,11,k+1)=0;endend%end%wucha=abs(u-uu);%求绝对误差矩阵11x11x141wucha1=wucha(:,:,11);%计算t=0.1时刻的绝对误差矩阵11x11wucha2=wucha(:,:,51

8、);%计算t=0.5时刻的绝对误差矩阵11x11wucha3=wucha(:,:,101);%计算t=1.0时刻的绝对误差矩阵11x11wucha4=wucha(:.:.141);%计算t=1.4时亥1I的绝对误差矢口阵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（'绝对误差值'）；%作t=0.5时刻的绝对误差图subpl

9、ot（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

10、,4）;mesh（x0,y0,wucha4）;title（'t=1.4时亥1J的绝对误差'）;xlabel（'x变量'）;ylabel（'y变量'）;zlabel（'绝对误差值'）；%1%1解%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

11、,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）;mesh（x0,y0,u（:,:,51）;%作t=0.5时刻的数值解title（'t=0.5时刻的数值解'）;xlabel（'x变量'）;ylabel（'y变量'）;zlabel（'u值'）；subplot（2,2

12、,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变量'）;ylabel（'y变量'）;zlabel（'

13、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（

14、'u值'）；subplot（2,2,4）;mesh（x0,y0,uu（:,:,141）;%作t=1.4时刻的精确解title（'t=1.4时刻的精确解'）;xlabel（'x变量'）;ylabel（'y变量'）;zlabel（'u值'）；五.实验结果及实例分析1、t=0.1、051.R1.4时刻的数值解与精确解图图1t=0.1、0.5时刻的数值解、精确解图2t=1.0、1.4时刻的数值解、精确解注：上两图为四个时刻的数值解与精确解，r='=0.1<-L（p代表维数）,本文p=2h、.P,三层显格式达二

15、阶收敛，不难看出，收敛效果很好，符合理论。下图是四个时刻的绝对误差图像，从图中看出，绝对误差较小，且经过计算得到，收敛阶近似于2,正好符合理论值2、t=0.1、051.01.4时刻的绝对误差图t=1。划的拾"一蝴弱的施好总XMUJM的推时-一mof”屁星灯盘.无0.02作。1序.0.1?10006.0.-3,1口图3四个时刻的绝对误差3、四个时刻（t=0.1、0.5、1.0、1.4）的绝对误差表t=0.1时刻的绝对误差0.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00010.0001

16、0.00020.00020.00020.00020.00020.00010.00010.00000.00000.00010.00030.00040.00040.00050.00040.00040.00030.00010.00000.00000.00020.00040.00050.00060.00060.00060.00050.00040.00020.00000.00000.00020.00040.00060.00070.00070.00070.00060.00040.00020.00000.00000.00020.00050.00060.00070.00080.00070.00060.0005

17、0.00020.00000.00000.00020.00040.00060.00070.00070.00070.00060.00040.00020.00000.00000.00020.00040.00050.00060.00060.00060.00050.00040.00020.00000.00000.00010.00030.00040.00040.00050.00040.00040.00030.00010.00000.00000.00010.00010.00020.00020.00020.00020.00020.00010.00010.00000.00000.00000.00000.0000

18、0.00000.00000.00000.00000.00000.00000.0000t=0.5时刻的绝对误差0.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00070.00130.00180.00210.00220.00210.00180.00130.00070.00000.00000.00130.00250.00340.00400.00420.00400.00340.00250.00130.00000.00000.00180.00340.00470.00550.00580.00550.0047

19、0.00340.00180.00000.00000.00210.00400.00550.00650.00680.00650.00550.00400.00210.00000.00000.00220.00420.00580.00680.00710.00680.00580.00420.00220.00000.00000.00210.00400.00550.00650.00680.00650.00550.00400.00210.00000.00000.00180.00340.00470.00550.00580.00550.00470.00340.00180.00000.00000.00130.0025

20、0.00340.00400.00420.00400.00340.00250.00130.00000.00000.00070.00130.00180.00210.00220.00210.00180.00130.00070.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000t=1.0时刻的绝对误差0.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00160.00310.00430.00510.00530.0051

21、0.00430.00310.00160.00000.00000.00310.00590.00820.00960.01010.00960.00820.00590.00310.00000.00000.00430.00820.01130.01320.01390.01320.01130.00820.00430.00000.00000.00510.00960.01320.01560.01640.01560.01320.00960.00510.00000.00000.00530.01010.01390.01640.01720.01640.01390.01010.00530.00000.00000.0051

22、0.00960.01320.01560.01640.01560.01320.00960.00510.00000.00000.00430.00820.01130.01320.01390.01320.01130.00820.00430.00000.00000.00310.00590.00820.00960.01010.00960.00820.00590.00310.00000.00000.00160.00310.00430.00510.00530.00510.00430.00310.00160.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000t=1.4时刻的绝