第三章線性方程組求解的數(shù)值方法

第三章線性方程組求解的數(shù)值方法3.1Gauss消去法與矩陣的LU分解3.2Cholesky分解3.3向量范數(shù)與矩陣范數(shù)3.4古典迭代法的構(gòu)造3.5迭代法的分析3.6超松弛迭代（SOR）及分塊迭代方法3.7線性方程組的條件3.8稀疏矩陣的計(jì)算線性方程組求解的數(shù)值方法§3.1Gauss消去法與矩陣的LU分解基本思想：用逐次消去未知數(shù)的方法把原方程組化為三角形方程組再求解。消元：用初等變換將原方程組的系數(shù)矩陣化為三角形矩陣（簡(jiǎn)稱(chēng)三角陣）再求解的方法。回代：解出三角形方程組的最后一個(gè)方程，將求得的值逐步往前一個(gè)方程代入的方法。消元為什么選主元避免方法：高斯主元消元法functionX=uptrbk(A,b)%A是一個(gè)n階矩陣。%b是一個(gè)n維向量。%X是線性方程組AX=b的解。[NN]=size(A);X=zeros(1,N+1);Aug=[Ab];forp=1:N-1[Y,j]=max(abs(Aug(p:N,p)));C=Aug(p,:);

Aug(p,:)=Aug(j+p-1,:);Aug(j+p-1,:)=C;

高斯主元消去法的MATLAB實(shí)現(xiàn)

ifAug(p,p)==0'A是奇異陣，方程無(wú)惟一解'breakendfork=p+1:Nm=Aug(k,p)/Aug(p,p);Aug(k,p:N+1)=Aug(k,p:N+1)-m*Aug(p,p:N+1);endend%這里用到程序8-1中定義的函數(shù)backsub來(lái)進(jìn)行回代。X=backsub(Aug(1:N,1:N)，Aug(1:N,N+1)；

A=[0.00123;-13.7124.623;-21.0725.643];>>b=[123]';>>uptrbk(A,b)比較計(jì)算結(jié)果：主元消元法-0.4904-0.05100.3675

用x=A\b-0.4904-0.05100.3675

一般消元法-0.4-0.099800.4functionX=backsub(A,b)%A是一個(gè)n階上三角非奇異陣。%b是一個(gè)n維向量。%X是線性方程組AX=b的解。n=length(b);X=zeros(n,1);X(n)=b(n)/A(n,n);fork=n-1:-1:1X(k)=(b(k)-A(k,k+1:n)*X(k+1:n))/A(k,k);End用回代法求解上三角線性方程組AX=B，其中A為非奇異。矩陣的LU分解functionX=lufact(A,b)%A為n階矩陣。%b是n維向量。%X是所求的AX=b的解。[N,N]=size(A);X=zeros(N,1);Y=zeros(N,1);C=zeros(1,N);R=1:N;forp=1:N-1[max1,j]=max(abs(A(p:N,p)));C=A(p,:);A(p,:)=A(j+p-1,:);A(j+p-1,:)=C;d=R(p);R(p)=R(j+p-1);R(j+p-1)=d;ifA(p,p)==0'A是奇異陣，方程組無(wú)惟一解'

breakendfork=p+1:N

mult=A(k,p)/A(p,p);

A(k,p)=mult;A(k,p+1:N)=A(k,p+1:N)-mult*A(p,p+1:N);endendY(1)=b(R(1));fork=2:N

Y(k)=b(R(k))-A(k,1:k-1)*Y(1:k-1);endX(N)=Y(N)/A(N,N);fork=N-1:-1:1

X(k)=(Y(k)-A(k,k+1:N)*X(k+1:N))/A(k,k);end§3.2Cholesky分解（平方根法)§3.3向量范數(shù)與矩陣范數(shù)向量的收斂矩陣范數(shù)常用的矩陣范數(shù)矩陣的收斂譜半徑§3.4三種迭代法雅可比迭代法的Matlab程序functionX=jacobi(A,b,P,delta,max1)%A是n維非奇異陣。%b是n維向量。%P是初值。%delta是誤差界。%max1是給定的迭代最高次數(shù)。%X為所求的方程組AX=b的近似解。N=length(b);fork=1:max1forj=1:N

X(j)=(b(j)-A(j,[1:j-1,j+1:N])*P([1:j-1,j+1:N]))/A(j,j);enderr=abs(norm(X'-P));P=X'；

if(err<delta)breakendendX=X';k,err;functionX=gseid(A,b,P,delta,max1)%A是n維非奇異陣。%b是n維向量。%P是初值。%delta是誤差界。%max1是給定的迭代最高次數(shù)。%X為所求的方程組AX=b的近似解。N=length(b);fork=1:max1forj=1:Nifj==1X(1)=(b(1)-A(1,2:N)*P(2:N))/A(1,1);

elseifj==NX(N)=(b(N)-A(N,1:N-1)*(X(1:N-1))')/A(N,N);else

X(j)=(b(j)-A(j,1:j-1)*X(1:j-1)-A(j,j