双曲方程基于某matlab的数值解法.doc
word双曲型方程基于MATLAB的数值解法一：一阶双曲型微分方程的初边值问题 准确解为 二：数值解法思想和步骤2.1:网格剖分为了用差分方法求解上述问题，将求解区域作剖分。将空间区间作等分，将时间区间作等分，并记。分别称和为空间和时间步长。用两簇平行直线将分割成矩形网格。2.2：差分格式的建立2.2.1：Lax-Friedrichs方法对时间、空间采用中心差分使得如此由上式得到Lax-Friedrichs格式截断误差为所以Lax-Friedrichs格式的截断误差的阶式令：如此可得差分格式为其传播因子为： 化简可得：所以当时，,格式稳定。* 2.2.2：LaxWendroff方法用牛顿二次插值公式可以得到LaxWendroff的差分格式，在此不详细分析，它的截断误差为，是二阶精度；当时，格式稳定。在这里主要用它与上面一阶精度的Lax-Friedrichs 方法进展简单比照。差分格式的求解因为时格式稳定，不妨取差分格式写成如下矩阵形式：如此需要通过对k时间层进展矩阵作用求出k+1时间层。对上面的矩阵形式通过matlab编出如附录的程序求出数值解、真实解和误差。 算法以与结果function P U E x t=PDEHyperbolic(uX,uT,M,N,C,type)format long%一阶双曲型方程的差分格式 %P U E x t=PDEHyperbolic(uX,uT,M,N,C,phi,psi1,psi2,type) %方程:u_t+C*u_x=0 0 <= t <= uT, 0 <= x <= uX%初值条件:u(x,0)=phi(x)%输出参数:U -解矩阵% x -横坐标% t -纵坐标，时间 % 输入参数:uX -变量x的上界% uT -变量t的上界% M -变量x的等分区间数% N -变量t的等分区间数% C -系数% phi -初值条件函数，定义为联函数% psi1,psi2 -边值条件函数，定义为联函数 % type -差分格式，从如下值中选取 % -type='LaxFriedrichs'，采用Lax-Friedrichs差分格式求解 % -type='LaxWendroff'，采用Lax-Wendroff差分格式求解 h=uX/M;%变量x的步长k=uT/N;%变量t的步长 r=k/h;%步长比 x=(0:M)*h; t=(0:N)*k; U=zeros(M+1,N+1); %初值条件for i=1:M+1 U(i,1)=cos(pi*x(i); P(i,1)=cos(pi*x(i); E(i,1)=0;end%边值条件for j=1:N+1 U(1,j)=cos(pi*t(j); E(1,j)=0; P(1,j)=cos(pi*t(j); U(M+1,j)=-cos(pi*t(j); P(M+1,j)=-cos(pi*t(j); E(M+1,j)=0;endswitch typecase'LaxFriedrichs'if abs(C*r)>1 disp('|C*r|>1，Lax-Friedrichs差分格式不稳定！') end%逐层求解 for j=1:N for i=2:M U(i,j+1)=(U(i+1,j)+U(i-1,j)/2-C*r*(U(i+1,j)-U(i-1,j)/2; P(i,j+1)=cos(pi*(x(i)+t(j+1); E(i,j+1)=abs(U(i,j+1)-cos(pi*(x(i)+t(j+1);endend%Lax-Wendroff差分格式 case'LaxWendroff'if abs(C*r)>1 disp('|C*r|>1，Lax-Wendroff差分格式不稳定！') end%逐层求解 for j=1:N for i=2:M U(i,j+1)=U(i,j)-C*r*(U(i+1,j)-U(i-1,j)/2+C2*r2*(U(i+1,j)-2*U(i,j)+U(i-1,j)/2; P(i,j+1)=cos(pi*(x(i)+t(j+1); E(i,j+1)=abs(U(i,j+1)-cos(pi*(x(i)+t(j+1);endendotherwise disp('差分格式类型输入有误！') return; endU=U'P=P'E=E'%作出图形 准确解mesh(x,t,P); title('一阶双曲型方程的准确解图像'); xlabel('空间变量 x'); ylabel('时间变量 t'); zlabel('一阶双曲型方程的解 P')%作出图形 数值解mesh(x,t,U); title(type '格式求解一阶双曲型方程的解的图像'); xlabel('空间变量 x'); ylabel('时间变量 t'); zlabel('一阶双曲型方程的解 U')return;命令窗口输入：>>uX=1;uT=1;M=90;N=100;C=-1;phi=inline('cos(pi*x)');psi1=inline('cos(pi*t)');psi2=inline('-cos(pi*t)');type='LaxFriedrichs'或type='LaxWendroff'>>P U E x t=PDEHyperbolic(uX,uT,M,N,C,type)从 matlab的数值解法结果中抽出一局部数据进展比拟表1LaxFriedrichs格式jkx,t数值解真实解误差4611(0.5,0.1)4621(0.5,0.2)4631(0.5,0.3)4641(0.5,0.4)4651(0.5,0.5)4661(0.5,0.6)4671(0.5,0.7)4681(0.5,0.8)4691(0.5,0.9)46101(0.5,1.0)表2LaxWendroff格式jkx,t数值解真实解误差4611(0.5,0.1)124621(0.5,0.2)-0204631(0.5,0.3)-0224641(0.5,0.4)0164651(0.5,0.5)-0014661(0.5,0.6)-0174671(0.5,0.7)0344681(0.5,0.8)0484691(0.5,0.9)05746101(0.5,1.0)6备注：本来，但是由于matlab中下标必须从大于0开始，所以在程序中图像分析:结果分析：从表1和表2可以看出LaxFriedrichs格式和LaxWendroff格式的真值得误差都比拟小，而LaxWendroff格式虽然精度比LaxFriedrichs的精度高，但是在网格点划分比拟细的情况下，二者的差异不大。从三个图像的结果看出，二者都拟合的相当好，并且结果都稳定。9 / 9
