《矩阵的LU分解(自编MATLAB)实验报告31893.pdf》由会员分享，可在线阅读，更多相关《矩阵的LU分解(自编MATLAB)实验报告31893.pdf（7页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、.1/7 1 矩阵的LU分解 1.1 LU 分解原理 定理：设ACnn，如果A的顺序主子式 0,0,，0 则存在唯一的主对角线上元素全为 1 的下三角矩阵L与唯一的上三角矩阵U，使得 A=LU.证明：对矩阵A的阶数使用数学归纳法.显然，当 n=1 时，=1就是唯一的分解式。现假定对 n-1 阶矩阵，定理的结论成立。对A进行分块 A=其中.由于 n-1 阶矩阵的k阶顺序主子式就是 A 的k阶主子式（k=1,2，n-2），故它们都不为零.从而由归纳法假设，有唯一的LU分解 其中的主对角线上的元素都 1.由于=0 所以是 n-1 阶可逆矩阵 先假设已有A=LU，其中 L=，U=.2/7，=A 则，是
2、 n-1 阶可逆矩阵，则由上式可惟一确定,这就证明了A的 LU分解的存在性和唯一性.1.2 LU 分解算法 当n 阶矩阵满足定理的条件时，可以用初等变换的方法求出L和 U.因为当A=LU时，由于L 可逆，故必存在可逆矩阵P 使得 即PA=PLU=U.也就是说，可以先对A施行行的初等变换得出上三角矩阵U，而矩阵P可以通过对单位矩阵I进行相同的行初等变换得出，即 P(A,I)(PA,PI)(U,P)于是,为保持P为下三角矩阵（从而也是下三角矩阵），在进行行初等变换时，不能进行行的对换，上行的倍数应加到下行的对应元.1.3 LU 分解用于解方程组 矩阵的三角分解在求解线性方程组时十分方便.如对线性方
3、程组 由于 是下三角矩阵，则解向量 可以通过依次求出其分量而求出，在求解方程组可以通过该方程组依次求出分量而快速得出.于是由两个方程组的求解而给出 的解.1.4 程序流程图 .3/7 输入矩阵A判断A是否为nn矩阵？计算A的n-1阶顺序主子式是否为0？输出：“无法进行LU分解”设定 n 阶单位矩阵L和全零矩阵 U输出L,U是是否否计算U 的第一行 U1j=A1j j=1，2，n计算L 的第一列U 的第 r 行（逐行算出）L 的第 r 列（逐列算出）1.5 MATLAB 程序 function f=LU_decom(A)m,n=size(A)if m=n.4/7 fprintf(Error:m
4、and n must be equal!m=%d,n=%dn,m,n)end for i=1:n-1 if(det(A(1:i,1:i)=0)fprintf(Error:det A(%d,%d)=0!n,i,i)flag=failure return;else flag=ok;end end L=eye(n);U=zeros(n);for i=1:n U(1,i)=A(1,i);end for r=2:n L(r,1)=A(r,1)/U(1,1);end for i=2:n for j=i:n z=0;for r=1:i-1 z=z+L(i,r)*U(r,j);end U(i,j)=A(i,j
5、)-z;end if abs(U(i,i)A=2 1 1;4 1 0;-2 2 1;LU_decom(A)m=3 n=3 L=1 0 0 2 1 0 -1 -3 1 U=2 1 1.6/7 0 -1 -2 0 0 -4（2）解方程组，程序与结果如下%-用LU分解解线性方程组-y=zeros(n,1);y(1)=b(1);for i=2:n y(i)=b(i)-sum(L(i,1:i-1).*y(1:i-1);end y x(n)=y(n)/U(n,n);for i=n-1:-1:1 x(i)=(y(i)-sum(U(i,i+1:n).*x(i+1)/U(i,i);end x=x 运行结果如下：y=1 0 2 x=-0.5000 1.0000 -0.5000 1.7 数据分析 调用MATLAB固有的LU分解函数，以与解方程组相关函数对以上数据进行计算，运行结果如下：A=2 1 1;4 1 0;-2 2 1;b=1 2 1;L,U=lu(A)L=0.5000 0.2000 1.0000.7/7 1.0000 0 0 -0.5000 1.0000 0 U=4.0000 1.0000 0 0 2.5000 1.0000 0 0 0.8000 X=inv(A)*b X=0.2500 1.0000 -0.5000 经比对结果相同，可见以上程序可行。