数值逼近上机练习答案.docx

《数值逼近上机练习答案.docx》由会员上传分享，免费在线阅读，更多相关内容在教育资源-天天文库。

1、实验一1．分别用循环型Arnoldi算法和循环型GMRES算法求解线性方程组。取，初始值，迭代终止条件为。要求输出数值近似解和迭代步数。M-文件Arnoldi算法functionX=Arnoldi(X0,m,A,b,ep)nk=0;while(max(abs(A*X0-b))>=ep)r0=b-A*X0;n=length(b);v=zeros(n,m+1);h=zeros(m,m);p=sqrt(sum(r0.^2));v(:,1)=r0./p;fork=2:m;tmp=0;fori=1:k-1h(i,k-1)=(A*v(:,k-1))'*v(:,i);tmp=

2、tmp+h(i,k-1)*v(:,i);endvk=A*v(:,k-1)-tmp;h(k,k-1)=sqrt(sum(vk.^2));v(:,k)=vk/h(k,k-1);endtmp=0;fori=1:mh(i,m)=(A*v(:,m))'*v(:,i);tmp=tmp+h(i,m)*v(:,i);endvm1=A*v(:,m)-tmp;hm1=sqrt(sum(vm1.^2));v(:,m+1)=vm1/hm1;e1=zeros(m,1);e1(1)=1;y=inv(h)*p*e1;z=v(:,1:m)*y;X0=X0+z;nk=nk+1;endnkX=X0

3、;Commandwindows：A=[1000;1100;1110;1111]b=[1111]'m=2x0=[0000]'ep=1.0e-8结果：Arnoldi(X0,m,A,b,ep)nk=23ans=1.00000.0000-0.0000-0.0000GMRES算法functionX=GMRES(X0,m,A,b,ep)nk=0;while(max(abs(A*X0-b))>=ep)r0=b-A*X0;n=length(b);v=zeros(n,m+1);h=zeros(m,m);p=sqrt(sum(r0.^2));v(:,1)=r0./p;fork=2:

4、m;tmp=0;fori=1:k-1h(i,k-1)=(A*v(:,k-1))'*v(:,i);tmp=tmp+h(i,k-1)*v(:,i);endvk=A*v(:,k-1)-tmp;h(k,k-1)=sqrt(sum(vk.^2));v(:,k)=vk/h(k,k-1);endtmp=0;fori=1:mh(i,m)=(A*v(:,m))'*v(:,i);tmp=tmp+h(i,m)*v(:,i);endvm1=A*v(:,m)-tmp;hm1=sqrt(sum(vm1.^2));v(:,m+1)=vm1/hm1;H=[h;[zeros(1,m-1),hm1

5、]];e1=zeros(m+1,1);e1(1)=1;y=pinv(H)*p*e1;z=v(:,1:m)*y;X0=X0+z;nk=nk+1;endnkX=X0;CommandwindowsA=[1000;1100;1110;1111]b=[1111]'m=2x0=[0000]'ep=1.0e-8结果：GMRES(X0,m,A,b,ep)nk=12ans=1.00000.0000-0.0000-0.00001.用追赶法、线性插值法和双参数法求解阶三对角方程组。其中，要求分别输出当时的计算误差和计算时间。M文件双参数法function[X,deta]=sdjscs

6、(a,b,c,n,f)if(length(f)>=n&length(f)<=n)s=zeros(n,1);t=zeros(n,1);tics(1)=0;s(2)=f(1)/c(1);t(1)=1;t(2)=-b(1)/c(1);fori=3:ns(i)=(f(i-1)-b(i-1)*s(i-1)-a(i-1)*s(i-2))/c(i-1);t(i)=-(b(i-1)*t(i-1)+a(i-1)*t(i-2))/c(i-1);endx=zeros(n,1);x(1)=(f(n)-a(n)*s(n-1)-b(n)*s(n))/(a(n)*t(n-1)+b(n)*t(

7、n));fori=2:nx(i)=t(i)*x(1)+s(i);endX=x;A=[zeros(1,n);[diag(a(2:end))zeros(n-1,1)]]+diag(b)+[[zeros(n-1,1)diag(c(1:end-1))];zeros(1,n)];deta=norm(f-A*X,2);tocelseX='youarewrong!';end窗口：n=1000;a=ones(1,n-1);b=-4*ones(1,n);c=2*ones(1,n-1);A=diag(b)+diag(a,-1)+diag(c,1);f=ones(n,1);f(1)=

8、3;f(end)=-3;