矩陣?yán)碚摲治鰉atlab實驗

1、矩陣?yán)碚搶嶒瀸I(yè)：控制科學(xué)與工程姓名：劉思學(xué)號：2115080100201 矩陣的LU分解1.1 原理定義1.1設(shè)ACn×n，若A可以表示成一個下三角矩陣L與一個上三角矩陣U的乘積A=LU（1.1）則稱其為矩陣A 的LU分解（三角分解）.定理1.1設(shè)ACn×n，如果A的順序主子式a110，a11a12a21a220，a11a1 n-1an-1 1an-1 n-10（1.2）則存在唯一的主對角線上元素全為一的下三角矩陣L與唯一的上三角矩陣U，使得A=LU. 曾祥金.矩陣分析及應(yīng)用.武漢大學(xué)出版社.2007.8：74-751.2 算法此算法較為簡單，直接敘述如下：由公式aij=

2、k=1nlikukj（1.3），可以按照從U的第一行求起，再求L的第一列，再求U的下一行，再求L的下一列，循環(huán)下去，直到求出LU.1.3 流程圖1.4 程序1 使用MATLAB的自帶程序進行LU分解A=2 1 1;4 1 0;-2 2 1;L,U,P = lu(A);運行結(jié)果A=2 1 1;4 1 0;-2 2 1;>> L,U,P=lu(A)L = 1.0000 0 0 -0.5000 1.0000 0 0.5000 0.2000 1.0000U = 4.0000 1.0000 0 0 2.5000 1.0000 0 0 0.8000P = 0 1 0 0 0 1 1 0 02

3、自編程序進行LU分解function L,U = LSLU( A )%UNTITLED 劉思編寫的矩陣的LU分解程序n,n=size(A); %獲取矩陣A的階L=zeros(n,n); %將矩陣LU置零U=zeros(n,n); for i=1:n %對L矩陣的對角線元素賦值1L(i,i)=1; end for k=1:n %根據(jù)算法求L,U矩陣for j=k:n U(k,j)=A(k,j)-sum(L(k,1:k-1).*U(1:k-1,j)'); end for i=k+1:n L(i,k)=(A(i,k)-sum(L(i,1:k-1).*U(1:k-1,k)')/U(k,

4、k); endendend運行結(jié)果：>> A=2 1 1;4 1 0;-2 2 1A = 2 1 1 4 1 0 -2 2 1>> A=2 1 1;4 1 0;-2 2 1A = 2 1 1 4 1 0 -2 2 1>> L,U=LSLU(A)L = 1 0 0 2 1 0 -1 -3 1U = 2 1 1 0 -1 -2 0 0 -4與MATLAB 自帶分解分解程序一致。2 矩陣的QR分解2.1 原理定理2.1設(shè)ACm×n，mn且rank A=n，則比存在非奇異上三角n×n矩陣R及m×n矩陣Q， QHQ=In，使得A=QR（2

5、.1）則稱其為矩陣A 的QR分解. 曾祥金.矩陣分析及應(yīng)用.武漢大學(xué)出版社.2007.8：78-792.2 算法1 利用 Schmidt正交化求矩陣的 QR 分解，Schmidt正交化方法是矩陣的 QR 分解最常用的方法. 主要依據(jù)下面的兩個結(jié)論 :結(jié)論 1 設(shè) A 是 n階實非奇異矩陣 ,則存在正交矩陣 Q 和實非奇異上三角矩陣 R 使 A 有 QR 分解 ;且除去相差一個對角元素的絕對值（模）全等于 1的對角矩陣因子外 ,分解是唯一的. 2設(shè) A 是 m × n實矩陣 ,且其n個列向量線性無關(guān) ,則 A 有分解 A =QR,其中 Q 是 m ×n實矩陣 ,且滿足 Q H

6、TQ = E,R 是n階實非奇異上三角矩陣該分解除去相差一個對角元素的絕對值 （模） 全等于1的對角矩陣因子外是唯一的。 李建東.矩陣 QR分解的三種方法.呂梁高等?？茖W(xué)校學(xué)報.2009.3:16步驟 : 1、寫出矩陣的列向量；2、把列向量按照 Schmidt正交化方法進行正交； 3、得出矩陣的 Q'R'；4、對Q'的列向量單位化得到Q，并在R'的每一行乘以Q'每一列模得到R.2.3 流程圖2.4 程序1.MATLAB自帶函數(shù)進行QR分解B=1 1;0 1;1 1;Q,R=qr(B)Q,R,E=qr(B)運行結(jié)果B=1 1;0 1;1 1;Q,R=qr(

7、B)Q,R,E=qr(B)Q = -0.7071 0 -0.7071 0 1.0000 0 -0.7071 0 0.7071R = -1.4142 -1.4142 0 1.0000 0 0Q = -0.5774 0.4082 -0.7071 -0.5774 -0.8165 0 -0.5774 0.4082 0.7071R = -1.7321 -1.1547 0 0.8165 0 0E =0 11 0自編算法進行QR分解function Q,R = LSQR( A )%LSQR 劉思編寫的矩陣的QR分解程序A=1,1;0,1;1,1;m,n=size(A); %得到A的維數(shù)if (m>=n

8、)&&(n=rank(A) R=eye(n);%R=E Q=zeros(m,n);%對Q賦初值 Q(:,1)=A(:,1); %b1=a1 for i=2:n %提取A的第i列，即A的第i個列向量 for j=1:i-1 x=Q(:,j); % bj R(j,i)=dot(A(:,i),x)/dot(x,x); %確立系數(shù)陣R end Q(:,i)=A(:,i)-Q(:,1:i-1)*R(1:i-1,i); end for i=1:n y=Q(:,i); R(i,:)= R(i,:)* norm(y); %對Q單位化（列向量單位化）； Q(:,i)= Q(:,i)/ norm(

9、y); %對R每一行乘以相應(yīng)的常數(shù)使得A=QR endelse disp('不滿足QR分解要求！')endend運行結(jié)果A = 1 1 0 1 1 1>> QQ = 1.0000 0 0 0 -0.7071 -0.7071 0 -0.7071 0.7071>> RR = 1.0000 1.0000 -0.7071 -1.4142 0.7071 0.0000與MATLAB 自帶分解分解程序一致。3 矩陣的奇異值分解3.1 原理設(shè)ACm×n，s1，s2，sr 是A 的非零奇異值，則存在m階酉矩陣U及n階酉矩陣V，m×n矩陣D.D= = 使

10、得A=UDVH這就是矩陣A的奇異值分解.3.2 算法第一步：求出AHA的特征值0=，確定非零奇異值=，i=1,2，r.第二步：分別求出矩陣AHA的對應(yīng)于特征值的特征向量并將其單位正交化，得到標(biāo)準(zhǔn)正交向量組1，2，n.令V=（1，2，n）=（V1，V2），V1=（1，2，r），V2=（r+1，r+2，m）.第三步：若U=（1，2，r，r+1，r+2，m）=（U1，U2）其中U1=（1，2，r），U2=（r+1，r+2，m），則因（A1，A2，Ar）=（A1，A2，Ar）即有U1=AV1 .其中=.第四步：解方程組AAHy=0，對基礎(chǔ)解系單位正交化可以求得r+1,r+2,m，令U=（1，2，r，r

11、+1，r+2，m）.曾祥金.矩陣分析及應(yīng)用.武漢大學(xué)出版社.2007.8：86-883.3 流程圖3.4 程序1. MATLAB自帶函數(shù)進行奇異值分解A=1 0;0 1;1 0;U,S,V=svd(A)%矩陣的奇異值分解運行結(jié)果A=1 0;0 1;1 0;U,S,V=svd(A)%矩陣的奇異值分解U = -0.7071 0 -0.7071 0 1.0000 0 -0.7071 0 0.7071S = 1.4142 0 0 1.0000 0 0V = -1 0 0 12自編程序進行奇異值分解function U,S,V = LSSVD(A)%LSQR 劉思編寫的矩陣的奇異值分解程序m,n=size(A); %得到A的維數(shù)U=zeros(m);V=zeros(n);r=rank(A); S=zeros(m,n);B,C=eig(A'*A);x=diag(C);B=B.',x;B=sortrows(B,-(n+1);for i=1:r S(i,i)=sqrt(B(i,n+1);endB=B(:,1:n);B=B.'V=qr(B);V1=V(:,1:r);U(:,1:r)=A*V1*(inv(S(1:r,1:r);U(:,r+1:m)=