微分方程数值解大作业

微分方程数值解大作业双曲型方程求解如下双曲型方程其中关于时间和空间的二阶偏导数均采用二阶中心差商逼近.1．，且初边条件如下：存在精确解为：2．，且初边条件如下：存在精确解为：第一题数值解的源程序为：functionu=hyperbolic(minx,maxx,mint,maxt,M,N,a)h=(maxx-minx)/M;%横轴步长t=(maxt-mint)/N;%纵轴步长r=t/h;u=zeros(M+1,N+1);%先定义u为M+1*N+1维的矩阵u(1,:)=0;%左边界为0u(M,:)=0;%右边界为0u(:,1)=bottoma(h,M)';%下边界由已知函数bottoma计算fori=1:M-1%对下边界上的导数g1的差分格式进行定义g1(i)=g(i*h);endforj=2:M%由导数的差分格式求第2行的u值u(j,2)=0.5*(a^2*r^2*(u(j+1,1)+u(j-1,1))+2*(1-a^2*r^2)*u(j,1)+t*g1(j-1));endfork=2:N%从第三行开始，一行一行地求u值forj=2:Mu(j,k+1)=a^2*r^2*(u(j+1,k)+u(j-1,k))+2*(1-a^2*r^2)*u(j,k)-u(j,k-1);endend下边界的导数的定义functiong=g(x)g=0下边界函数的定义functiona=bottoma(h,m)pi=3.1415926;a(1)=0;a(m+1)=0;fori=2:ma(i)=sin(pi*(i-1)*h);end精确解的源程序：functionu=hy(minx,maxx,mint,maxt,M,N)pi=3.1415926;x=minx+(maxx-minx)*(0:M)/M;t=mint+(maxt-mint)*(0:N)/N;forj=1:N+1fori=1:M+1u(i,j)=cos(pi*t(j))*sin(pi*x(i));endend稳定时的误差sum(sum(hyperbolic(0,1,0,2,4,8,1)-hy(0,1,0,2,4,8)))ans=2.587544328334521e-007误差很小第二题数值解代码functionu=hyperbolic(minx,maxx,mint,maxt,M,N,a)h=(maxx-minx)/M;%横轴步长t=(maxt-mint)/N;%纵轴步长r=t/h;u=zeros(M+1,N+1);%先定义u为M+1*N+1维的矩阵u(1,:)=0;%左边界为0u(M,:)=0;%右边界为0u(:,1)=bottoma(h,M)';%下边界由已知函数bottoma计算fori=1:M-1g1(i)=g(i*h);%对下边界上的导数g1的差分格式进行定义endforj=2:M%由导数的差分格式求第2行的u值u(j,2)=0.5*(a^2*r^2*(u(j+1,1)+u(j-1,1))+2*(1-a^2*r^2)*u(j,1)+t*g1(j-1));endfork=2:N%从第三行开始，一行一行地求u值forj=2:Mu(j,k+1)=a^2*r^2*(u(j+1,k)+u(j-1,k))+2*(1-a^2*r^2)*u(j,k)-u(j,k-1);endend下边界函数的定义functiona=bottoma(h,m)a=0;下边界的导数的定义functiong=g(x)pi=3.141592653589793;g=sin(4*pi*x);稳定的精确解functionu=hy(minx,maxx,mint,maxt,M,N)pi=3.141592653589793;x=minx+(maxx-minx)*(0:M)/M;t=mint+(maxt-mint)*(0:N)/N;forj=1:N+1fori=1:M+1u(i,j)=sin(t(j))*sin(4*pi*x(i));endend稳定性条件和Courant条件分析。我们采用分离变量的方法对格式进行稳定性分析，可得传播矩阵：G而这个矩阵的特征方程为：λ他的根λ1,2按模不大于1的充要条件为2-ar≤1.亦即，当上式成立时，vonNeuma