矩阵LU分解求逆详细分析与C语言实现

1、题目要求给定一个多维矩阵，实现该矩阵的求逆运算。1、理论分析矩阵的一种有效而广泛应用的分解方法是矩阵的LU三角分解，将一个n阶矩阵A分解为一个下三角矩阵L和一个上三角矩阵U的乘积。所以首先对矩阵进行三角分解，这里采用Doolittle分解，即分解为一个下三角矩阵（对角元素为1），和一个上三角矩阵的乘积。再进行相应的处理。所以，矩阵求逆的算法流程可表述如下：图1 矩阵求逆流程图1)进行LU分解；2)对分解后的L阵（下三角矩阵）和U阵（上三角矩阵）进行求逆;；3)L阵的逆矩阵和U阵的逆矩阵相乘，即可求得原来矩阵的逆。即： （1）1.1矩阵的LU 分解若n阶方阵 的各阶顺序主子式不等于零，即： （2

2、）则A的LU分解存在且唯一。 （3）由矩阵的乘法原理, 可推导出LU分解的迭代算法（4） （5） （6） （7）矩阵的LU分解是一个循环迭代的过程, U矩阵是从第1行迭代到第n行, 而L矩阵则是从第1列迭代到第n列, 且U矩阵先于L矩阵一个节拍。1.2 L矩阵和U矩阵求逆首先假设下三角矩阵L的逆矩阵为，不失一般性，考虑4阶的情况，利用，有：（1） ，；（2）（3）（4）。从而求得下三角矩阵L的逆矩阵R式如下：， （8）上三角矩阵U的逆矩阵可以由下式得到：。， （9）矩阵求逆是一个迭代的过程，依次循环, 迭代次, 求出整个逆矩阵。其中U矩阵的循环迭代时按行顺序，列倒序进行，L矩阵的循环迭代按列顺

3、序，行顺序进行，直到计算出整个矩阵的所有结果为止。1.3 矩阵相乘上三角矩阵U的逆矩阵u与下三角矩阵L的逆矩阵相乘, 最终得到原始矩阵A的逆矩阵, 完成整个矩阵求逆的过程。对于n阶矩阵相乘的迭代形式可表示如下： （10）1.4 实例分析例：给定一4阶矩阵，通过LU分解求逆矩阵。解：算法过程为：，第一步：求LU矩阵设，通过（4）（7）式可逐步进行矩阵L和U中元素的计算，如下所示：经迭代计算，最后得到L和U矩阵为：第二步：求L和U矩阵的逆,（1）求U矩阵的逆由式（9）可得矩阵U的逆的各元素计算如下：（2）求L矩阵的逆由（8）式可得L矩阵的逆的各元素计算如下所以得到L和U的逆矩阵为：（3）求A的逆矩

4、阵由式（10）可计算得到矩阵A的逆，如下：由程序计算出的结果如下：2、C语言程序设计及测试2.1 算法c程序实现#include#include #define N 4void main() float aNN; float LNN,UNN,outNN, out1NN; float rNN,uNN; memset( a , 0 , sizeof(a); memset( L , 0 , sizeof(L); memset( U , 0 , sizeof(U); memset( r , 0 , sizeof(r); memset( u , 0 , sizeof(u); int n=N; int k

5、,i,j; int flag=1; float s,t; /input a matrix/ printf(ninput A=); for(i=0;in;i+) for(j=0;jn;j+) scanf(%f,&aij);/figure the input matrix/printf(输入矩阵:n);for(i=0;in;i+)for (j = 0; j n; j+)printf(%lf , aij);printf(n);for(j=0;jn;j+) a0j=a0j; /计算U矩阵的第一行 for(i=1;in;i+) ai0=ai0/a00; /计算L矩阵的第1列 for(k=1;kn;k+)

6、 for(j=k;jn;j+) s=0; for (i=0;ik;i+) s=s+aki*aij; /累加 akj=akj-s; /计算U矩阵的其他元素 for(i=k+1;in;i+) t=0; for(j=0;jk;j+) t=t+aij*ajk; /累加 aik=(aik-t)/akk; /计算L矩阵的其他元素 for(i=0;in;i+) for(j=0;jj) Lij=aij; Uij=0;/如果ij，说明行大于列，计算矩阵的下三角部分，得出L的值，U的/为0 else Uij=aij; if(i=j) Lij=1; /否则如果ij，说明行小于列，计算矩阵的上三角部分，得出U的/值，

7、L的为0 else Lij=0; if(U11*U22*U33*U44=0) flag=0; printf(n逆矩阵不存在); if(flag=1) /求L和U矩阵的逆 for (i=0;i=0;k-)s=0;for (j=k+1;j=i;j+)s=s+Ukj*uji;uki=-s/Ukk;/迭代计算，按列倒序依次得到每一个值，for (i=0;in;i+) /求矩阵L的逆 rii=1; /对角元素的值，直接取倒数，这里为1for (k=i+1;kn;k+)for (j=i;j=k-1;j+)rki=rki-Lkj*rji; /迭代计算，按列顺序依次得到每一个值/绘制矩阵LU分解后的L和U矩阵

8、/ printf(nLU分解后L矩阵:); for(i=0;in;i+) printf(n); for(j=0;jn;j+) printf( %lf,Lij); printf(nLU分解后U矩阵:); for(i=0;in;i+) printf(n); for(j=0;jn;j+) printf( %lf,Uij); printf(n); /绘制L和U矩阵的逆矩阵 printf(nL矩阵的逆矩阵:); for(i=0;in;i+) printf(n); for(j=0;jn;j+) printf( %lf,rij); printf(nU矩阵的逆矩阵:); for(i=0;in;i+) prin

9、tf(n); for(j=0;jn;j+) printf( %lf,uij); printf(n); /验证将L和U相乘，得到原矩阵 printf(nL矩阵和U矩阵乘积n); for(i=0;in;i+) for(j=0;jn;j+)outij=0; for(i=0;in;i+) for(j=0;jn;j+) for(k=0;kn;k+) outij+=Lik*Ukj; for(i=0;in;i+) for(j=0;jn;j+) printf(%lft,outij); printf(rn); /将r和u相乘，得到逆矩阵 printf(n原矩阵的逆矩阵:n); for(i=0;in;i+) for(j