1. 用matlab編寫一個Gauss消去法的程序求解方程組=程序：>> clear>> A=2 1 0;3 4 -1;5 2 6A = 2 1 0 3 4 -1 5 2 6>> b=1;0;1b = 1 0 1>> L,U=lu(A)L = 0.4000 0.0714 1.0000 0.6000 1.0000 0 1.0000 0 0U = 5.0000 2.0000 6.0000 0 2.8000 -4.6000 0 0 -2.0714>> c=Lbc = 1.0000 -0.6000 0.6429>> x=Ucx = 0.8621 -0.7241 -0.31032、列主元Gauss消去法求解方程組=程序：A=input(請輸入線性方程組的增廣矩陣A=);n=length(A)-1;x=zeros(n,1);aa=zeros(n,1);for j=1:nfor i=1:(n+1)AA(j,i)=abs(A(j,i);endendfor k=1:(n-1)for i=k:naa(i-(k-1)=AA(i,k);endfor i=k:nif AA(i,k)=max(aa)breakendendif AA(i,k)=0breakfprintf(方程組系數(shù)矩陣奇異n);elsefor j=k:(n+1)jh=A(i,j);A(i,j)=A(k,j);A(k,j)=jh;endendfenzi=A(k,k);for j=k:(n+1)A(k,j)=A(k,j)/fenzi;endfor p=(k+1):njj=A(p,k);for j=k:(n+1)A(p,j)=A(p,j)-jj*A(k,j);endendendif k=(n-1)x(n)=A(n,(n+1)/A(n,n);for i=(n-1):(-1):1he=0;for j=(i+1):nhe=he+A(i,j)*x(j);endx(i)=A(i,(n+1)-he;endendx請輸入線性方程組的增廣矩陣A=2 4 0 1;3 -5 1 1;1 7 8 1x = 0.4032 0.04840.03233、用LU分解求解線形方程組=程序：>> A=2 1 -1 5;1 6 4 -1;3 2 7 5;0 -1 8 2A = 2 1 -1 5 1 6 4 -1 3 2 7 5 0 -1 8 2>> b=8;4;2;1b = 8 4 2 1>> L,U=lu(A)L = 0.6667 -0.0625 -0.6692 1.0000 0.3333 1.0000 0 0 1.0000 0 0 0 0 -0.1875 1.0000 0U = 3.0000 2.0000 7.0000 5.0000 0 5.3333 1.6667 -2.6667 0 0 8.3125 1.5000 0 0 0 2.5038>> L,U,P=lu(A)L = 1.0000 0 0 0 0.3333 1.0000 0 0 0 -0.1875 1.0000 0 0.6667 -0.0625 -0.6692 1.0000U = 3.0000 2.0000 7.0000 5.0000 0 5.3333 1.6667 -2.6667 0 0 8.3125 1.5000 0 0 0 2.5038P = 0 0 1 0 0 1 0 0 0 0 0 1 1 0 0 0>> y=L(P*b)y = 2.0000 3.3333 1.6250 7.9624>> x=Uyx = -5.3063 2.3333 -0.37843.18024、用Cholesky分解求解線形方程組=程序：>> clear>> A=6 3 -2;3 5 1;-2 1 7A = 6 3 -2 3 5 1 -2 1 7>> b=2;-1;1b = 2 -1 1>> L=chol(A)L = 2.4495 1.2247 -0.8165 0 1.8708 1.0690 0 0 2.2783>> Y=LbY = 0.8165 -1.0690 1.2332>> X=LYX = 0.9541 -0.8807 0.54135、利用Jacobic迭代法求解方程組=，eps=1.010，=程序：function x, k, index=Jacobi(A, b, eps, it_max)if nargin <4 it_max=100; endif nargin <3 eps; endn=length(A); k=0; x=zeros(n,1); y=zeros(n,1); index=1;while 1 for i=1:n y(i)=b(i); for j=1:n if j=i y(i)=y(i)-A(i,j)*x(j); end end if abs(A(i,i)<1e-10 | k=it_max index=0; return; end y(i)=y(i)/A(i,i); end if norm(y-x,inf)<eps break; end x=y; k=k+1;end運行：>> A=0.98 -0.05 -0.02;-0.04 -0.9 0.07;-0.02 0.09 0.94；>> b=1;1;1；>> x,k,index=Jacobi(A,b,1e-6,100)；x = 0.9904 -1.0628 1.1867k = 6index = 16、利用G-S迭代法求解方程組=，esp=1.010,=程序：function y=gaussseidel(A,b,x0,esp)D=diag(diag(A);U=-triu(A,1);L=-tril(A,-1);G=(D-L)U;f=(D-L)b;y=G*x0+f;n=1;while norm(y-x0)>=esp&n<=1000 x0=y; y=G*x0+f; n=n+1;endyn答案：A=10 3 1;2 -10 3;1 3 10; b=14;-5;14; x0=0;0;0; esp=1.0e-3; gaussseidel(A,b,x0,esp)y = 1.00003896860354 1.00002773079501 0.99998778390114n = 6ans = 1.00003896860354 1.00002773079501 0.999987783901147、利用SOR迭代法求解方程組 =Eps=1.0*10-6,取程序：function n,x=sor(A,b,X,nm,w,ww)%用超松弛迭代法求解方程組Ax=b%輸入：A為方程組的系數(shù)矩陣，b為方程組右端的列向量，X為迭代初值構(gòu)成的列向量，nm為最大迭代次數(shù)，w為誤差精度，ww為松弛因子%輸出：x為求得的方程組的解構(gòu)成的列向量，n為迭代次數(shù)n=1;m=length(A);D=diag(diag(A); %令A(yù)=D-L-U,計算矩陣DL=tril(-A)+D; %令A(yù)=D-L-U,計算矩陣LU=triu(-A)+D; %令A(yù)=D-L-U,計算矩陣UM=inv(D-ww*L)*(1-ww)*D+ww*U); %計算迭代矩陣g=ww*inv(D-ww*L)*b; %計算迭代格式中的常數(shù)項%下面是迭代過程while n<=nm x=M*X+g; %用迭代格式進行迭代 if norm(x-X,inf)<w disp(迭代次數(shù)為);n disp(方程組的解為);x return; %上面：達到精度要求就結(jié)束程序，輸出迭代次數(shù)和方程組的解 end X=x;n=n+1;end%下面：如果達到最大迭代次數(shù)仍不收斂，輸出警告語句及迭代的最終結(jié)果（并不是方程組的解）disp(在最大迭代次數(shù)內(nèi)不收斂!);disp(最大迭代次數(shù)后的結(jié)果為);x運行：>> A=0.78 -0.02 -0.12 -0.14;-0.02 0.86 -0.04 0.06;-0.12 -0.04 0.72 -0.08;-0.14 0.06 -0.08 0.74; b=0.76;0.08;1.12;0.68;c=1000; d=1e-6; f=1.03; k=0;0;0;0;g=sor(A,b,k,c,d,f)迭代次數(shù)為n = 7方程組的解為x = 1.5350 0.1220 1.9752 1.4130g = 7 8、用共軛梯度法求解上題的線性方程組，Eps=1.0*10-6,9、利用冪法求：的按模最大的特征值和對應(yīng)的特征向量，限定最大的迭代步驟n=500，Eps=1.0，=程序：function k,lambda,Vk,Wc=mifa(A,x0,eps,n)lambda=0;k=1;Wc =1; state=1; V=x0;while(k<=n)&(state=1)Vk=A*V; m j=max(abs(Vk); mk=m; tzw=abs(lambda-mk); Vk=(1/mk)*Vk;Txw=norm(V-Vk); Wc=max(Txw,tzw); V=Vk;lambda=mk;state=0;if(Wc>eps)state=1;endk=k+1;Wc=Wc;endif(Wc<=eps)disp(迭代次數(shù)k,主特征值的近似值lambda,主特征向量的近似向量Vk,相鄰兩次迭代的誤差Wc如下：) elsedisp(迭代次數(shù)k已經(jīng)達到最大迭代次數(shù)max1,主特征值的迭代值lambda,主特征向量的迭代向量Vk,相鄰兩次迭代的誤差Wc如下：) end Vk=V;k=k-1;Wc;主程序：>> A=1 2 3;2 3 4;3 4 5;x0=1;1;1;esp=1.0e-6;n=500;k,lambda,xk,Wc=mifa(A,x0,eps,n),k,lambda,xk,Wc=mifa(A,x0,0.000001,500),x,D = eig (A), Dzd=max(diag(D), wuD= abs(Dzd- lambda), wux=x(:,2)./xk,答案：迭代次數(shù)k已經(jīng)達到最大迭代次數(shù)max1,主特征值的迭代值lambda,主特征向量的迭代向量Vk,相鄰兩次迭代的誤差Wc如下：k = 500lambda = 9.6235xk = 0.5247 0.7623 1.0000Wc = 1.7764e-015迭代次數(shù)k,主特征值的近似值lambda,主特征向量的近似向量Vk,相鄰兩次迭代的誤差Wc如下：k = 8lambda = 9.6235xk = 0.5247 0.7623 1.0000Wc = 1.5190e-007x = 0.8277 0.4082 0.3851 0.1424 -0.8165 0.5595 -0.5428 0.4082 0.7339D = -0.6235 0 0 0 -0.0000 0 0 0 9.6235Dzd = 9.6235wuD = 9.2421e-009wux = 0.7781 -1.0710 0.408210、利用反冪法求A=的按模最小的特征值和對應(yīng)的特征向量，限定最大的迭代步驟n=500，Eps=1.0，=程序：function k,lambdan,Vk,Wc=ydwyfmf(A,V0,jlamb,jd,max1)n,n=size(A); A1=A-jlamb*eye(n); jd= jd*0.1;RA1=det(A1); if RA1=0disp(請注意：因為A-aE的n階行列式hl等于零，所以A-aE不能進行LU分解.)returnendlambda=0;if RA1=0 for p=1:nh(p)=det(A1(1:p, 1:p);endhl=h(1:n);for i=1:nif h(1,i)=0disp(請注意：因為A-aE的r階主子式等于零，所以A-aE不能進行LU分解.) returnendend if h(1,i)=0 disp(請注意：因為A-aE的各階主子式都不等于零，所以A-aE能進行LU分解.)k=1;Wc =1;state=1; Vk=V0;while(k<=max1)&(state=1)L U=lu(A1); Yk=LVk;Vk=UYk; m j=max(abs(Vk); mk=m;Vk1=Vk/mk; Yk1=LVk1;Vk1=UYk1; m j=max(abs(Vk1); mk1=m;Vk2=(1/mk1)*Vk1;tzw1=abs(mk-mk1)/mk1);tzw2=abs(mk1-mk);Txw1=norm(Vk)-norm(Vk1);Txw2=(norm(Vk)-norm(Vk1)/norm(Vk1);Txw=min(Txw1,Txw2); tzw=min(tzw1,tzw2); Vk=Vk2;mk=mk1; Wc=max(Txw,tzw); Vk=Vk2;mk=mk1;state=0;if(Wc>jd)state=1;endk=k+1;%Vk=Vk2,mk=mk1,endif(Wc<=jd)disp(A-aE的秩R(A-aE)和各階順序主子式值hl、迭代次數(shù)k,按模最小特征值的近似值lambda,特征向量的近似向量Vk,相鄰兩次迭代的誤差Wc如下：) elsedisp(A-aE的秩R(A-aE)和各階順序主子式值hl、迭代次數(shù)k已經(jīng)達到最大迭代次數(shù)max1,按模最小特征值的迭代值lambda,特征向量的迭代向量Vk,相鄰兩次迭代的誤差Wc如下：) endhl,RA1endendV,D=eig(A,nobalance),Vk;k=k-1;Wc;lambdan=jlamb+1/mk1;主程序：A=1 2 0;4 -3 2;5 1 4;x0=1;1;1;k,lambda,Vk,Wc=ydwyfmf(A,x0,0.2,0.000001,500)請注意：因為A-aE的各階主子式都不等于零，所以A-aE能進行LU分解.A-aE的秩R(A-aE)和各階順序主子式值hl、迭代次數(shù)k,按模最小特征值的近似值lambda,特征向量的近似向量Vk,相鄰兩次迭代的誤差Wc如下：hl = 0.8000 -10.5600 -21.7280RA1 = -21.7280V = -0.1577 -0.5531 0.3776 -0.3243 -0.0508 -1.0000 -1.0000 1.0000 -0.1070D = 5.1127 0 0 0 1.1836 0 0 0 -4.2964k = 7lambda = 1.1836Vk = 0.5531 0.0508 -1.0000Wc = 5.3678e-009