《偏微分实验报告三》由会员分享，可在线阅读，更多相关《偏微分实验报告三（10页珍藏版）》请在金锄头文库上搜索。
1、重 庆 大 学学 生 实 验 报 告实验课程名称 偏微分方程数值解 开课实验室 数统学院 学 院 数统 年级 2012 专业班 信计1班 学 生 姓 名 张伟 学 号 20122058 开 课 时 间 2014 至 2015 学年第 2 学期总 成 绩教师签名数学与统计学院制推荐精选开课学院、实验室： 数统学院 实验时间 ： 2015年 5月20日实验项目名 称两点边值问题的有限差分法实验项目类型验证演示综合设计其他指导教师李茂军成 绩是一实验目的通过该实验，要求学生掌握求解两点问题的有限差分法，并能通过计算机语言编程实现。二实验内容考虑如下的初值问题： (1) (2)其中，是给定常数。 将区
2、间等分，设，网点。 1在第三部分写出问题（1）和（2）的差分格式，并给出该格式的局部截断误差。2根据你写出的差分格式，编写一个有限差分法程序。将所写程序放到第四部分。3给定参数，问题(1)的精确解，其中将及带入方程(1)可得。分别取，用所编写的程序计算问题(1)和(2)。将数值解记为，网点处精确解记为，。然后计算相应的误差推荐精选，及收敛阶，将计算结果填入第五部分的表格，并对表格中的结果进行解释？4. 将数值解和精确解画图显示，每种网格上的解画在一张图。要求：将程序放到实验结果部分。将电子版的实验报告发送到邮箱,邮件标题为：偏微分方程实验报告2-姓名-学号。纸质版的实验报告在下次上机实验时提交
3、。三实验原理、方法（算法）、步骤用差商代替导数的方法：设为内网点,用表示方括号内的函数在出的取值。由Taylor展开有： （3）类似地， （4）为简化符号，记，由（3）和（4），以及边值问题的解在网点满足，故 （5）其中与的导数有关的量 （6）若舍去，则得到逼近边值问题（1）和（2）的差分格式为 （7）推荐精选 （8）其中称为差分算子，称为差分方程（7）的局部截断误差，利用，将（5）写为 （9） 将（7）改写为： （10）其中 （11）这是以为未知量的线性代数方程。由于，当充分小时，必有。另外，当，有。若记方程组的未知向量,右端向量，系数矩阵为 （12） 则有 （13） 易知，为不可约对角占优
4、矩阵，故为非奇异的，即上述线性方程组有唯一解。由于为三对角矩阵，故用追赶法求解。四实验环境（所用软件、硬件等）及实验数据文件Matlab1、两点边值问题的有限差分法程序clear;a=0;b=1;%定义初始值alpha=0;beta=0;%定义初始u值N=161;%这里的N代表一共多少个网点（包括边界点），题中要求，10,20,40,80,160，对应数字加1推荐精选h=(b-a)/(N-1); %步长x=linspace(0,1,N);N=N-1; %二阶线性常微分方程两点边值问题系数确定for i=1:N-1 r(i)=1;%u(x)一阶导系数r(x)的确定 p(i)=3;%变系数p(x)
5、的确定 q(i)=2;%u(x)系数函数q（x） f(i)=(-10*x(i+1)-1)*exp(x(i+1);%f(x)的确定endfor i=1:N-1 b1(i)=2*(p(i)+p(i)/h+2*h*q(i);%b（i）,如果上式p(i)不是常数，需重新定义endfor i=1:N-2 c1(i)=2*p(i)/h-r(i); %c（i），p（i+1/2）endfor i=1:N-2 a1(i)=2*p(i)/h+r(i); %a（i）,p（i-1/2）end%写为HU=G的形式H=diag(b1)+diag(-a1,-1)+diag(-c1,1);for i=1:N-1 g1(i)=
6、2*h*f(i);%这里考虑到alpha=0;belta=0，如不是，需如下定义end %g1(1)=2*h*f(1)+a1(1)*alpha; %g1(N-1)=2*h*f(N-1)+c1(N-1)*beta;G=g1; %右端矩阵GY=inv(H)*G;U=zeros(1,N+1);U(1)=alpha;U(N+1)=beta;推荐精选for i=1:N-1 U(i+1)=Y(i);end plot(x,U,r-*)title(数值解与精确解对比图);xlabel(x值);ylabel(u值);legend()grid onhold onR=x.*(x-1).*exp(x); plot(x
7、,R,g+-);wucha=abs(R-U);hold onplot(x,wucha,b*-)legend(精确解u,数值解u,绝对误差u);eN_C=max(wucha)eN_0=sqrt(h*sum(wucha.2)2、收敛阶的计算clc;eN_C=0.0013 3.3514e-04 8.3970e-05 2.0994e-05 5.2486e-06;eN_0=9.7257e-04 2.4314e-04 6.0867e-05 1.5218e-05 3.8045e-06;for i=1:4 SLJC(i)=(log(eN_C(i)/eN_C(i+1)/log(2); SLJ0(i)=(log(
8、eN_0(i)/eN_0(i+1)/log(2); fprintf( 收敛阶SLJC(%d)分别为： %fn,i,SLJC(i) fprintf( 收敛阶SLJ0(%d)分别为： %fnn,i,SLJ0(i)end推荐精选推荐精选五实验结果及实例分析1、收敛阶计算结果N收敛阶收敛阶10 0.0013-9.7257e-04-203.3514e-041.9556762.4314e-042.000015408.3970e-051.9968186.0867e-051.998055802.0994e-051.9998971.5218e-051.9998811605.2486e-061.9999733.8045e-062.000000由表中数据不难看出，收敛阶趋于二阶，即符合。2、数值解与精确解比较图如下： 图1 h=0.025数值解与精确解比较注：其它N=10,20,80,160,的图像都大致相同，这里不再累述，图片已保存在压缩包内。 教师签名年 月 日 （注：可编辑下载，若有不当之处，请指正，谢谢!） 推荐精选
