资源描述：
《六点对称格式.doc》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。
1、利用差分法求解下列问题：方向，方向上。在观察数值解与精确解的误差。一、算法描述1.网格剖分取，方向，方向上。2.差分格式选用六点对称格式，即Crank-Nicolson格式，是将向前差分格式和向后差分格式做算术平均得到的一种差分格式。，，，，得到矩阵形式：3.初边值处理直接给出。在本问题中，边值都为0，所以程序简化了很多，编程时和都直接取为0。一、程序%*****************清理缓存*******************************clearall;clc;dx=0.1;%空间步长
2、dt=0.01;%时间步长t0=0;tn=0.25;%时间起点；时间终点x0=0;xn=1;%空间起点；空间终点N=(xn-x0)/dx;%空间网格数T=(tn-t0)/dt;%时间网格数fori=1:N+1x(i)=(i-1)*dx;end%****************边界条件*********************************%u(0,t)=0,u(1,t)=0;%可知在边界上的值始终为0，后面的时间每层迭代都直接取为0.%****************初始条件**********
3、***********************u=sin(pi*x)+x.*(1-x);r=dt/dx^2;%****************储存三对角元素***************************%%为减小储存量，只储存三对角元素%A=zeros(3,N-1);fori=1:N-1A(1,i)=-r/2;A(1,N-1)=0;A(2,i)=1+r;A(3,1)=0;A(3,i)=-r/2;endA=mylu(A);%LU分解a=r/2;b=1-r;c=r/2;forj=1:T%*******
4、***计算方程组右边的值**************d(1)=b*u(2)+c*u(3)+0.02;fori=2:N-2d(i)=a*u(i)+b*u(i+1)+c*u(i+2)+0.02;endd(N-1)=a*u(N-1)+b*u(N)+0.02;%d=myzhuigan(A,d);%追赶法求解u=[0d0];%两个0是给出的边界条件endplot(x,u,'o')%绘制数值解图形holdu_xt=exp(-pi*pi*T*dt).*sin(pi.*x)+x.*(1-x);plot(x,u_xt,'r
5、')%绘制精确解图形e=u_xt-u%误差%**************over****************************************子程序1：function[LU]=mylu(A)%这里的LU分解是为了分解只储存三对角元素的矩阵%%forexample%%A=[22220;21111;0-1-1-1-1];%%%result:%LU=%2.00002.00002.00002.00000%2.00002.00002.00002.00002.0000%0-0.5000-0.500
6、0-0.5000-0.5000%n=size(A,2);%矩阵的列数fori=2:nA(3,i)=A(3,i)/A(2,i-1);A(2,i)=A(2,i)-A(3,i)*A(1,i-1);endLU=A;子程序2：function[x]=myzhuigan(A,d)n=length(d);%*****************追的过程**********************************fori=2:nd(i)=d(i)-A(3,i)*d(i-1);end%*****************
7、*赶的过程*********************************d(n)=d(n)/A(2,n);fori=n-1:-1:1d(i)=(d(i)-A(1,i)*d(i+1))/A(2,i);endx=d;三、结果：
