数值逼近机练习答案.docx

实验一分别用循环型Arnoldi算法和循环型GMRES算法求解线性方程组。取，初始值，迭代终止条件为。要求输出数值近似解和迭代步数。M-文件Arnoldi算法function X=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;for k=2:m; tmp=0;for i=1:k-1 h(i,k-1)=(A*v(:,k-1))*v(:,i); tmp=tmp+h(i,k-1)*v(:,i);end vk=A*v(:,k-1)-tmp; h(k,k-1)=sqrt(sum(vk.^2)); v(:,k)=vk/h(k,k-1);endtmp=0;for i=1:m h(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; Command windows：A=[1 0 0 0;1 1 0 0;1 1 1 0;1 1 1 1]b=[1 1 1 1]m=2x0=[0 0 0 0]ep=1.0e-8结果：Arnoldi(X0,m,A,b,ep)nk = 23ans = 1.0000 0.0000 -0.0000 -0.0000GMRES算法function X=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;for k=2:m; tmp=0;for i=1:k-1 h(i,k-1)=(A*v(:,k-1))*v(:,i); tmp=tmp+h(i,k-1)*v(:,i);end vk=A*v(:,k-1)-tmp; h(k,k-1)=sqrt(sum(vk.^2)); v(:,k)=vk/h(k,k-1);endtmp=0;for i=1:m h(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]];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; Command windowsA=[1 0 0 0;1 1 0 0;1 1 1 0;1 1 1 1]b=[1 1 1 1]m=2x0=[0 0 0 0]ep=1.0e-8结果：GMRES(X0,m,A,b,ep)nk = 12ans = 1.0000 0.0000 -0.0000 -0.0000用追赶法、线性插值法和双参数法求解阶三对角方程组。其中，要求分别输出当时的计算误差和计算时间。M文件双参数法function [X,deta]=sdjscs(a,b,c,n,f)if(length(f)=n length(f)=n) s=zeros(n,1); t=zeros(n,1); tic s(1)=0; s(2)=f(1)/c(1); t(1)=1; t(2)=-b(1)/c(1);for i=3:n s(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);end x=zeros(n,1); x(1)=(f(n)-a(n)*s(n-1)-b(n)*s(n))/(a(n)*t(n-1)+b(n)*t(n));for i=2:n x(i)=t(i)*x(1)+s(i);end X=x; A=[zeros(1,n);[diag(a(2:e