《矩阵LU分解求逆详细分析与C语言实现.docx》由会员分享,可在线阅读,更多相关《矩阵LU分解求逆详细分析与C语言实现.docx(7页珍藏版)》请在优知文库上搜索。
1、题目要求给定一个多维矩阵,实现该矩阵的求逆运算。1、理论分析矩阵的一种有效而广泛应用的分解方法是矩阵的1.U三角分解,将一个n阶矩阵A分解为一个下三角矩阵1.和一个上三角矩阵U的乘积。所以首先对矩阵进行三角分解,这里采用DOolittIe分解,即分解为一个下三角矩阵(对角元素为1),和一个上三角矩阵的乘积。再进行相应的处理。所以,矩阵求逆的算法流程可表述如下:图1矩阵求逆流程图1)进行1.U分解;2)对分解后的1.阵(下三角矩阵)和U阵(上三角矩阵)进行求逆;;3)1.阵的逆矩阵和U阵的逆矩阵相乘,即可求得原来矩阵的逆。即:A-,=(1.U)T=UT1.(1)1.1矩阵的1.U分解假设n阶方阵
2、的各阶顺序主子式不等于零,即:a1/12aka2a22a2kr.a=.0,(%=1,2,),(2)akak2akk那么A的1.U分解A=1.XU存在且唯一。由矩阵的乘法原理,可推导出1.U分解的迭代算法UOj=4)/,(/=0,1,2,.,h-1),IiO=i-0,l,2,11-l),wOOr-%-EIlikUkjZ-I%=%一*%,Irrk=l(r=0,1,2,H1;7=r.,H1),r-lairEjiMkr1 k=l(7)(r=0,1,2,1,111;i=r+1,/21)矩阵的1.U分解是一个循环迭代的过程,U矩阵是从第1行迭代到第n行,而1.矩阵那么是从第1列迭代到第n歹J,且U矩阵先于
3、1.矩阵一个节拍。1.2 1.矩阵和U矩阵求逆首先假设下三角矩阵1.的逆矩阵为I,不失一般性,考虑4阶的情况,利用=/,有:(1)100=4)1.l=耳,,22=5;,,33=4;4o=T0O(411.lo);(3),20=,00(,21乙10+,22心20);4Go=-a(4140+4240+公及)。从而求得下三角矩阵1.的逆矩阵R式如下:/Tjg7=1-(),j%5,ij0矩阵求逆是一个迭代的过程,依次循环,迭代-1次,求出整个逆矩阵。其中U矩阵的循环迭代时按行顺序,列倒序进行,1.矩阵的循环迭代按列顺序,行顺序进行,直到计算出整个矩阵的所有结果为止。1.3 矩阵相乘上三角矩阵U的逆矩阵U
4、与下三角矩阵1.的逆矩阵/相乘,最终得到原始矩阵A的逆矩阵AT=UT1.T=R,完成整个矩阵求逆的过程。对于n阶矩阵相乘的迭代形式可表示如下:(10)*皿刃团伙x伙刃k=j1.4实例分析4215通过1.U分解求逆矩阵A-例:给定一4阶矩阵A=:;?48366849解:算法过程为:A-=(1.XST=UX1.=IIX第一步:求1.U矩阵-1.oOo0。一UOOUOl002U031.XU=J1.IO1.Il0X0UnUnU3设乙201.?Ii22。00U22U23,通过(4)_乙304243_OOOU33(7)式可逐步进行矩阵1.和U中元素的计算,如下所示:经迭代计算,最后得到1.和U矩阵为:第二
5、步:求1.和U矩阵的逆,/(1)求U矩阵的逆由式(9)可得矩阵U的逆的各元素计算如下:(2)求1.矩阵的逆4k一.-000-1-1000-I,00000A1002100-I,10I11002201210V,20220由11.32%1.50.66671.251_,3041Iyl/33_式可得1.矩阵的逆的各元素计算如下所以得到1.和U的逆矩阵为:求A的逆矩阵由式(10)可计算得到矩阵A的逆,如下:A-1=ul0.25-0.166667-0.125-4.5100000.3333300-2100000.5-23-2100004_-1.9166670.83333-1.251.由程序8.833334-3
6、.666675.5-4.5-0.666670.33333005.333333-2.666673-2-7.6666673.33333-54计算出的结果如下:2、C语言程序设计及测试2.1算法C程序实现#include#include#defineN4voidmain()floataNN;float1.NN,UNN,outNN,outUNN;floatrNN,uNN;memset(a,0,sizeof(八));memset(1.,0,sizeof(1.);memset(U,0,sizeof(U);memset(r,0,sizeof(r);memset(u,0,sizeof(u);intn=N;in
7、tk,i,j;intflag=1;floats,t;/inputamatrix/printf(,ninputA=m);for(i=0;in;i+)for(j=0;jn;j+)scanf(,%f&aij);/figuretheinputmatrix/Printf(输入矩阵:n)for(i=0;in;i+)(for(j=0;jn;j+)(p11ntf(11%lf11,alij);)printf(11nn);)for(j=0;jn;j+)alj=alj;计算U矩阵的第一行for(i=l;in;i+)ai0=ai0a00;计算1.矩阵的第1列for(k=l;kn;k+)(for(j=k;jn;j+)(
8、s=0;for(i=0;ik;i+)s=s+aki*aig;累加ak11j=akj-s;计算U矩阵的其他元素for(i=k+l;in;i+)(t=0;for(j=0;jk;j+)t=t+aij*ajk;累加aik=(aik-t)akk;计算1.矩阵的其他元素)fdr(i=0;in;i+)for(j=0;jj)皿=0;如果ij,说明行大于列,计算矩阵的下三角局部,得出1.的值,U的为0elseUifj=aij;if(i=j)1.ij=l;否那么如果ivj,说明行小于列,计算矩阵的上三角局部,得出U的/值,1.的为0else1.ij=0;if(Ull*U22*U33*U44=0)(flag=0;P
9、rintf(n逆矩阵不存在);if(flag=l)/求1.和U矩阵的逆for(i=0;ivn;i+)/*求矩阵U的逆*/uii=lUii;对角元素的值,直接取倒数for(k=i-l;k=0;k-)s=0;for(j=k+hj=iy+)s=s+Ukj*uji;uk11i=sUkk;迭代计算,按列倒序依次得到每一个值,for(i=0;in;i+)/求矩阵1.的逆rii=l;对角元素的值,直接取倒数,这里为1for(k=i+l;kn;k+)for(j=i;j=k-l;j+)rki=rki-1.kj*rji;迭代计算,按列顺序依次得到每一个值/绘制矩阵1.U分解后的1.和U矩阵/PrindrW1.U分
10、解后1.矩阵:);for(i=0;in;i+)printf(11n);for(j=0;jn;j+)printf(,%l,1.ij);)PrindrW1.U分解后U矩阵:);for(i=0;in;i+)printf(,n11);for(j=0;jn;j+)P1111tf(11%lf,Uij);)printf(,n11);/绘制1.和U矩阵的逆矩阵Printfn1.矩阵的逆矩阵:);for(i=0;in;i+)printf(,n);for(j=0;jn;j+)printf(11%l,rij);)Printf(11U矩阵的逆矩阵巧;fbr(i=0;in;i+)printf(,n*);for(j=0;
11、jn;j+)P11ntf(,%lf,uij);Iprintf(11nu);验证将1.和U相乘,得到原矩阵Printf(n1.矩阵和U矩阵乘积n);for(i=0;in;i+)forG=0;jn;j+)outij=0;)for(i=0;in;i+)(for(j=0;jn;j+)(for(k=0;kn;k+)outifj=1.ik*Ukj;)for(i=0;in;i+)(for(j=0;jn;j+)(printf(%lt,outiirj);Prindrrn);)/将和U相乘,得到逆矩阵Printlrn原矩阵的逆矩阵:n)for(i=0;i、整数矩阵2、小数矩阵(2)满秩矩阵1整数矩阵的测试2小数矩阵的测试for(j=0;jn;j+)outlij=0;)for(i=0;in;i+)(for(j=0;jn;j+)(for(k=0;kn;k+)outlij+=uik*rkj;)for(i=0;in;i+)(for(j=0;jn;j+)(printf(11%lt11,outlirj);)printf(11rn11);)