上机矩阵实验讲述.docVIP

矩阵的LU分解 原理?： 设A若A可以表示成一个下三角矩阵和一个上三角矩阵的乘积A=LU，则称其为矩阵A的LU分解（三角分解）。? 矩阵的LU分解在求解线性方程组时将十分简便。如对线性方程组Ax=b，设A=LU是其LU分解。我们先求解方程组Ly=b。由于L是下三角矩阵，则解向量y可以通过依次求出其分量y1，y2，……，yn而求出，再求解方程组Ux=y。解向量x可以通过该方程组依次求出分量xn，xn-1，……，x2，x1而快速得出。于是由两个方程组Ux=y，Ly=b的求解而给出LUx=Ly=b=Ax的解。? 若矩阵A非奇异，则A能分解为LU的充分必要条件是A的顺序主子行列式不为0。?? ??? 则存在惟一的主对角线上元素全为1的下三角阵L与惟一的上三角阵U，使得A=LU。? 2、算法? 当n阶矩阵A满足LU分解的条件时，我们可以用初等变换的方法求出L和U.因为当A=LU时，由于L可逆，故必存在可逆阵P组PL=I 即 PA=PLU=U.也就是说，我们可以先对A施行行的初等变换得出上三角阵U,而矩阵P可以通过对单位阵I进行相同的行初等变换得出。即有P(A,I)=(PA,PI)=(U,P)于是A=P^-1U,为保持P为上三角阵，在进行行初等变换时，不能进行行的对换，上行的倍数应加到下行的对应元。 ? 我们知道对于第i行j列元素a[ij]的第k+1次变换，有 ? 设A=L*U，那么有如下数学表达： ? A=[;;;]; b=[1 1 1 1]; n=length(b);%方程个数n x=zeros(n,1);%未知向量 A(2:n,1)=A(2:n,1)./A(1,1); for i=2:n-1 A(i,i)=A(i,i)-sum(A(i,1:i-1).*A(1:i-1,i)); for j=i+1:n A(i,j)=A(i,j)-sum(A(i,1:i-1).*A(1:i-1,j)); A(j,i)=(A(j,i)-sum(A(j,1:i-1).*A(1:i-1,i)))/A(i,i); end end A(n,n)=A(n,n)-sum(A(n,1:n-1).*A(1:n-1,n)); A U=A;L=A; for i=1:n L(i,i)=1; end for i=1:n-1 for j=i+1:n L(i,j)=0; end end L %下三角阵 for i=2:n for j=1:i-1 U(i,j)=0; end end U %上三角阵 4、例子? 以矩阵A=[2?1?1；4?1?0；-2?2?1]为例说明上述程序的正确性，在matlab中输入以下指令：? ? 5、结果? 输入上述指令以后，得到如下的结果：? ? 此时可以说，上述程序是正确的。 矩阵分解成LU形式是有条件的，首先矩阵必须是非奇异的矩阵，其次矩阵的全部顺序主子式非零的时候才能完全保证矩阵可分解成LU且分解唯一。 矩阵的奇异值分解 1、原理 设A∈Cm×n，s1，s2，…，sr是A的非零奇异值，则存在m阶酉矩阵U∈Cm×n及n阶酉矩阵V，m×n矩阵D， D= = 使得A=UDVH这就是矩阵A的奇异值分解。 2、算法 第一步：求出AHA的特征值≥≥…≥＞0==…=，确定非零奇异值=，i=1,2，…，r。 第二步：分别求出矩阵AHA的对应于特征值的特征向量并将其单位正交化，得到标准正交向量组α1，α2，…，αn令V=（α1，α2，…，αn）=（V1，V2），V1=（α1，α2，…，αr），V2=（αr+1，αr+2，…，αn） 第三步：若U=（γ1，γ2，…，γr，γr+1，γr+2，…，γm）=（U1，U2），其中U1=（γ1，γ2，…，γr），U2=（γr+1，γr+2，…，γm）， 则因（Aα1，Aα2，…，Aαr）=（s1γ1，s2γ2，…，srγr） 即有U1=AV1 。其中= 第四步：解方程组AAHy=0，对基础解系单位正交化可以求得γr+1，γr+2，…，γm，令U=（γ1，γ2，…，γr，γr+1，γr+2，…，γm）。 3、程序 在matlab中有求解矩阵奇异值分解的svd函数，调用格式为[U，S，V] = svd(A)，其中U就是所求的U矩阵，S是所求的对角阵，V就是所求的酉矩阵V。 另一程序： function [U1,D,V1]=QYF(A)% 函数[L,D,V]=QYF(S)实现对矩阵A的奇异值分解 [m,n]=size(A); U1=zeros(m); V1=zeros(n); D=zeros(m,n); [V1,