资源描述：
《六点对称法,ADI法,预校法,和LOD法解二维抛物线方程.doc》由会员上传分享，免费在线阅读，更多相关内容在教育资源-天天文库。
1、偏微分数值解法实验报告实验名称：六点对称格式，ADI法，预校法，LOD法解二维抛物线方程的初值问题实验成员：吴兴杨敏姚荣华于潇龙余凡郑永亮实验日期：2013年5月17日指导老师：张建松一、实验内容用六点对称格式，ADI法，预校法和LOD法求解二维抛物线方程的初值问题：已知（精确解为：）设差分解为，则边值条件为：初值条件为：取空间步长，时间步长网比。1:ADI法：由第n层到第n+1层计算分成两步：先先第n层到n+1/2层，对uxx用向后差分逼近，对uyy用向前差分逼近，对uyy用向后差分逼近，于是得到了如下格式：其中j,k=1,2,…,M-1,n=0,1,2
2、,…,上标n+1/2表示在t=tn+1/2+(n+1/2)取值。假定第n层的已求得，则由（1）求出，再由（2）求出。2:预-校法差分格式：先通过U的n层求解U的n+1/4层，在通过U的n+1/4层求U的n的n+1/2层，最后通过U的n+1/2层求解U的n+1层，下为计算的预算格式：3:LOD算法：由第n层到第n+1层计算分为两步：（1）第一步：，构造出差分格式为：（2）第二步：，构造出差分格式为：其中。假定第n层的已求得，则由求出，这只需按行解一些具有三对角系数矩阵的方程组；再由求出，这只需按列解一些具有三对角系数矩阵的方程组，所以计算时容易实现的。4:六
3、点对称格式：将向前差分格式和向后差分格式做算术平均，即可以得到得六点对称格式：二、程序代码1:ADI法：%交替方向差分格式ADIclcx_a=0;x_b=1;%x的区间端点y_a=0;y_b=1;%y的区间端点N=40;%控制空间区域划分h=1/N;%空间步长x=[x_a:h:x_b];y=[y_a:h:y_b];T=1600;tao=1/T;%时间步长r=tao/(h^2);%网比a=1/16;U=ones(N+1,N+1);%迭代矩阵%按题意将边界点的值取为0forj=1:N+1U(1,j)=0;U(N+1,j)=0;end%初值条件fori=2:Nf
4、orj=1:N+1U(i,j)=sin(pi*x(i))*cos(pi*y(j));endend%差分格式方程组的系数矩阵diag_0=(1+r*a)*ones(N-1,1);diag_1=(-r*a/2)*ones(N-2,1)';A=diag(diag_0)+diag(diag_1,1)+diag(diag_1,-1);%组装系数矩阵A2=zeros(N+1);A2(2:N,2:N)=A;A2(1,1)=1;A2(N+1,N+1)=1;A2(1,2)=-1;A2(N+1,N)=-1;A2(2,1)=-r*a/2;A2(N,N+1)=-r*a/2;f=z
5、eros(N-1,1);f2=zeros(N+1,1);forn=1:T%计算到时间层t=1%x方向的迭代fork=2:Nforj=1:N-1%边界值为0，不必特殊处理j=1和N-1的情况f(j)=r*a/2*(U(j,k)+U(j+2,k))+(1-r*a)*U(j+1,k);endU(2:N,k)=Af;end%y方向的迭代forj=2:Nfork=2:Nf2(k)=r*a/2*(U(j,k+1)+U(j,k-1))+(1-r*a)*U(j,k);endU(j,:)=(A2f2)';endend%构造t=1时精确解网格函数jingquejie=ze
6、ros(N+1,N+1);fori=1:N+1forj=1:N+1jingquejie(i,j)=sin(pi*x(i))*cos(pi*y(j))*exp(-pi^2/8);endenddeta=abs(U-jingquejie);%绝对误差deta_max=max(max(deta));fprintf('最大误差%f',deta_max)figure(1);[x_l,y_l]=meshgrid(x);%生成网格采样点mesh(x_l,y_l,deta);title('误差网格分布');figure(2);mesh(x_l,y_l,jingqueji
7、e');%精确值的网格函数值title('精确解');figure(3);mesh(x_l,y_l,U');%数值解的网格函数title('数值解');U;2:预-校法差分格式：%用预-校法解抛物型方程clearclcformatlongJ=40;%x,y方向上的划分个数N=1600;%t方向上的划分个数，这里只求到t=1h=1/J;%x和y方向上的步长t=1/N;%t方向上的步长r=1;%网格比a=1/16;%方程中的系数[U]=zeros(J+1,J+1,N+1);%使用预-校法计算值[U1]=zeros(J+1,J+1,N+1);%真值%计算真值fo
8、rn=1:N+1fori=1:J+1forj=1:J+1U1(i,
