偏微分中心差分格式实验报告(含matlab程序).docx

二阶常微分方程的中心差分求解学校:中国石油大学（华东）理学院 姓名：张道德一、 实验目的1、 构造二阶常微分边值问题:其中为上的连续函数，为给定常数的中心差分格式；2、 根据中心差分格式求解出特定例题的数值解，并与该例题的精确解进行比较。二、 中心差分格式的构造将区间分成等分，分点为： 。于是我们得到区间的一个网络剖分。称为网格的节点称为步长。得中心差分格式为：其中式中。三、 差分格式求解根据中心差分格式可以构造出：可以看出系数矩阵为三对角矩阵，而对于系数矩阵为三对角矩阵的方程组可以用“追赶法”求解，则可以得出二阶常微分方程问题的数值解。四、 举例求解我们选取的二阶常微分方程边值问题为： 其精确解为：。则我们具体求解出的解如下：1、 当时，数值解与精确解为： (1) 表1、时，数值解与精确解统计表x的值0.10.20.30.40.5u的数值解1.0110691.0427441.0969041.1768961.28789u的精确解1.010051.0408111.0941741.1735111.284025两者之差0.0010190.0019340.0027290.0033850.003864x的值0.60.70.80.9u的数值解1.4374431.6363631.9000012.250209u的精确解1.4333291.6323161.8964812.247908两者之差0.0041140.0040460.003520.002301将两者绘于同一图中如下： (2)结论：可以看出数值解与精确解之间的误差很小，在 这样一个数量级上。我们也可以求出的范数 来，如下：Norm1()=0.0269;Norm2()=0.0095;Normoo()=0.0041;所以，可以得出中心差分格式求解该方程效果挺好。2、当时，将数值解与精确解绘于同一图像中，如下：(2)结论：可以看出数值解与精确解之间的误差很小，在 这样一个数量级上。我们也可以求出的范数 来，如下：Norm1()=0.0027;Norm2()=;Normoo()=;所以，可以得出中心差分格式求解该方程效果挺好。五、 程序程序1%*%f221.mfunction q,f=f211(x)%q函数，f函数q=4*x2;f=-2*exp(x2);%*程序2%*%追赶法function x=zhuiganfa(a,b,c,d)%对角线下方的元素，个数比A少一个% %对角线元素%对角线上方的元素，个数比A少一个%d为方程常数项%用追赶法解三对角矩阵方程r=size(a);m=r(2);r=size(b);n=r(2);if size(a)=size(c)|m=n-1|size(b)=size(d) error(变量不匹配,检查变量输入情况!);end%LU分解u(1)=b(1);for i=2:n l(i-1)=a(i-1)/u(i-1); u(i)=b(i)-l(i-1)*c(i-1); v(i-1)=(b(i)-u(i)/l(i-1); end%追赶法实现%求解Ly=d，追的过程y(1)=d(1);for i=2:n y(i)=d(i)-l(i-1)*y(i-1);end%求解Ux=y,赶的过程x(n)=y(n)/u(n);for i=n-1:-1:1 x(i)=y(i)/u(i); x(i)=(y(i)-c(i)*x(i+1)/u(i);end%*程序3%*%ODE2.mfunction x=ODE2(x0,xN,u0,uN,N)%中心差分求解%x0,xN初始条件%u0,uN边值条件%N等分数%步长h=1/N;%a(1:N-2)=-1/h2; %对角线下方的元素，个数比A少一个for i=1:N-1 z(i)=x0+i*h; q,f=f211(z(i); b(i)=2/h2+q; %对角线元素 d(i)=f; end%对对角线元素进行调整d(1)=d(1)+u0/h2;d(N-1)=d(N-1)+uN/h2;%x=zhuiganfa(a,b,a,d);%数值解%*程序4%*%main_chapter.m%x0,xN初始条件%u0,uN边值条件%N等分数%h步长x0=0;xN=1;u0=1;uN=exp(1);N=10;