大型矩阵快速求逆算法的研究.doc

《大型矩阵快速求逆算法的研究.doc》由会员分享，可在线阅读，更多相关《大型矩阵快速求逆算法的研究.doc（112页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、Four short words sum up what has lifted most successful individuals above the crowd: a little bit more.-author-date大型矩阵快速求逆算法的研究丽学院团2009 号说明书矩阵是高等代数的重要组成部分，是许多数学分支研究的重要工具。并且矩阵作为数学工具之一有其重要的使用价值，它常见于很多学科中，如：线性代数、线性规划、统计分析，以及组合数学等。在实际生活中，很多问题都可以借用矩阵抽象出来进行表述并进行运算，如在各循环赛中常用的赛况表格等，矩阵的概念和性质相对矩阵的运算较容易理解和掌握。

2、在矩阵的运算和应用中，仅逆矩阵的求解方法及算法问题就值得我们去好好研究，尤其是对大型矩阵的逆矩阵的研究。近年来，随着互联网的高速发展，计算机内部的运算量也急剧增加，如何把浩瀚的数据准确地计算出结果，并且加快它的计算速度，己成为一个备受关注的研究课题。随着计算机应用领域的发展，逆矩阵运算的需求越来越大。现阶段大型逆矩阵运算都是由软件实现，如Matlab等数学软件。逆矩阵运算软件的普及，将使计算机的逆矩阵运算性能得到几何级数的提高。伴随矩阵作、初等变换等算法作为逆矩阵运算的一个重要组成部分，对其研究的意义也就很大。 我们所做的研究方向，仅是利用Matlab数学软件和我们所学到的求逆矩阵的知识编写几

3、种求逆矩阵的算法。为了使内容更完整充实，文中列举了几种基本的求逆矩阵的方法以及整理了一些特殊矩阵求逆的一些简便方法等。最后简要说明了一般矩阵的广义逆矩阵。由于考虑自身学识的不足，我们并未对广义逆矩阵做详细的展开。希望有兴趣的读者能自己探索。大型矩阵快速求逆算法的研究信息09本耀魂雨大型矩阵快速求逆算法的研究本文首先介绍了逆矩阵的定义以及逆矩阵的相关性质。其次，根据逆矩阵的相关理论主要介绍了矩阵求逆的几种常用方法。 如定义法、伴随矩阵法、初等变换法、三角分解法、分块矩阵法等,并运用软件matlab7.0对一些方法实现了程序化。且通过多次检验证明了所编程序的正确性。文章最后简要阐述了对一些特殊矩阵

4、的求逆算法（如对称正定矩阵、有理矩阵），还有针对一般矩阵的广义逆矩阵的作了浅层次的探究。本文的相关研究不仅对提高矩阵求逆的速度和准确度起着较为重要的作用。而且对逆矩阵应用的进一步发展有着深远的意义。关键词：逆矩阵 求逆算法 Matlab7.0 广义矩阵逆矩阵的定义：对于n阶方阵A，如果有一个n 阶方阵B，使AB = BA = E，则说方阵 A 是可逆的，并把方阵 B 称为 A 的逆矩阵（其中只有方阵才有逆矩阵的概念）。矩阵可逆的条件： 定理1 若方阵 A 可逆，则 A 的行列式不等于 0 。 定理2 若A的行列式不等于0 ，则A可逆，且故，求逆矩阵可以用逆矩阵的定义来求或者利用求伴随矩阵法的方

5、法，还可以用行初等变换法（本文仅以行初等变换为例）、矩阵的三角分解、矩阵分块等方法求得逆矩阵。一、 五种基本求逆矩阵的方法1.1定义法求逆矩阵用一个例子来解释定义法求逆矩阵的方法：例：求矩阵A的逆矩阵，其中解：因为，所以存在.设=，由定义知,所以.由矩阵乘法得.由矩阵相等可解得故1.2伴随矩阵求逆矩阵（附算法）用求伴随矩阵的方法求逆矩阵的原理分析如下：设A=()为n阶矩阵，为|A|中元素的代数余子式, (i,j = 1, 2, ,n),则称矩阵=nnnnnnAAAAAAAAAALLLLLLLL212221212111*为A的伴随矩阵。根据上述定理2:若A的行列式不等于0,则A可逆，且, 1.3

6、初等变换法求逆矩阵（附算法）用初等变换法的方法求逆矩阵的原理分析如下：引理1.3.1 对m x n阶矩阵A，施行一次初等行变换，相当于在A的左边乘以相应m阶初等矩阵；对A施行一次初等列变换，相当于在A的右边乘以相应的n阶初等矩阵。公理1.3.1 初等矩阵都是可逆矩阵，其逆矩阵还是初等矩阵。推论1.3.1 A可逆的充要条件为A可表为若干初等矩阵之积。即推论1.3.2 A可逆，则A 可由初等行变换化为单位矩阵。 （1）由矩阵初等变换的这些性质可知，若A可逆，构造分块矩阵(AE)，其中E为与A同阶的单位矩阵，那么 （2）由（1）式 代入（2）式左边，有上式说明分块矩阵(AE)经过初等行变换，原来A的

7、位置变换为单位阵E，原来E的位置变换为我们所要求的,即1.4三角分解求逆矩阵（附算法）直接可以得到从此方法把A矩阵进行LU分解，并得到A的逆矩阵。1.5 分块矩阵法求逆此方法可将阶数较高的矩阵化为阶数较低的矩阵再求其逆，使计算简化1、 利用分块矩阵的乘法，求可逆矩阵的逆阵.设分块矩阵其中A，C均可逆，求解：令，于是有可得所以，从而2、 利用分块矩阵的初等变换，降阶求逆定理1.5.1 设A和B分别是可逆矩阵，C 是kr 矩阵，D 是rk 矩阵,“可逆”或“ 可逆”,则矩阵可逆。上例说明定理3可以应用到任何阶(奇数阶或偶数阶)可逆矩阵的求逆问题中。通常它使2n 阶可逆阵的求逆问题转化为一些n 阶可

8、逆阵的求逆问题，而使2m+1 阶可逆阵的求逆问题转化为一些m 阶和一些m+1 阶可逆阵的求逆问题，从而使计算的难度降低。因而这种求逆矩阵的方法也称之为“降阶法”。应当指出，当A,B 可逆时，相平行地，对于形如的矩阵(这里A,B 分别是k 阶和r 阶可逆矩阵,C 是kr 矩阵，D 是rk 矩阵）有着类似的结论：二、算法要求及算法分析：这里仅运用Matlab 7.0软件对求伴随矩阵、初等行变换及直接三角分解法求逆矩阵这三种方法编写了算法，并进行了适当的算法分析。2.1算法要求首先，不管是用求伴随矩阵A*还是利用初等变换求A的逆矩阵，都要先判定A矩阵是否为方阵。通过判断A矩阵的行数和列数是否相等来判

9、断A矩阵是否为方阵。其次，在确定A矩阵为方阵后还须判别A是否可逆。这里我们可以利用定理1，通过判别A的行列式是否为0，来确定A是否存在逆矩阵。最后，运用求伴随矩阵、初等行变换及直接三角分解法的计算方法来编写代码。2.2 三种求逆算法及解释分析以下是运用Matlab7.0编写的求逆矩阵的三个函数，以及对这三个函数的分析及比较：1、 利用求伴随矩阵的方法编写的求逆算法（代码如下）M函数 函数名：F=solvenif1(A)function F=solvenif1(A)n=rank(A);x,y=size(A);if x=y %判断A矩阵是否为方阵disp(请注意：因为该矩阵非方阵，所以该矩阵无逆.

10、)F=A无逆矩阵;return;else%disp(请注意：该矩阵为方阵.)d=det(A); if d=0 %判断A矩阵是否存在逆矩阵disp(请注意：因为该矩阵非满秩，所以该矩阵无逆.)F=A无逆矩阵;return; else %disp(请注意：因为该矩阵满秩，所以该矩阵逆存在.) for i=1:n for j=1:n B=A; B(i,:)=;B(:,j)=; C(j,i)=(-1)(i+j)*det(B)/d; end end F=C; disp(该矩阵的逆F= ) endend分析：利用求伴随矩阵的方法编写的求逆算法最主要的就是求元素的代数余子式。故一定要把元素所在的行和列删除，

11、且需要剩余部分组成的余子式输出。解释：B(i,:)=;B(:,j)=; %是为了选取元素的代数余子式。其中（i，j = 1,2,3n）为的代数余子式，B为（n-1）阶方阵。C(j,i)=(-1)(i+j)*det(B)/d; %是求A*与A的行列式d的比，结果即为（=/ d）。运行代码： clear;A=2 7 -1 9;-11 3 -3 7,F=solvenif1(A)A = 2 7 -1 9 -11 3 -3 7请注意：因为该矩阵非方阵,所以该矩阵无逆.F =A无逆矩阵 clear;A=2 3 -1;-1 3 -3;1 15 -11,F=solvenif1(A)A = 2 3 -1 -1

12、3 -3 1 15 -11请注意：因为该矩阵非满秩,所以该矩阵无逆.F =A无逆矩阵 clear;A=1 2 3;2 2 1;3 4 3;F=solvenif1(A)该矩阵的逆F= F = 1.0000 3.0000 -2.0000 -1.5000 -3.0000 2.5000 1.0000 1.0000 -1.0000 A*Fans = 1 0 0 0 1 0 0 0 1注：该算法对任何大型可逆方阵都适用2、 利用初等行变换求逆矩阵（代码如下）M函数 函数文件名：F=solvenif2(A)function F=solvenif2(A)n=rank(A);x,y=size(A);if x=y

13、 disp(请注意：因为该矩阵非方阵，所以该矩阵无逆.) F=A无逆矩阵;else %disp(请注意：该矩阵为方阵.) d=det(A); if d=0 disp(请注意：因为该矩阵非满秩，所以该矩阵无逆.) F=A无逆矩阵; else %disp(请注意：因为该矩阵满秩，所以该矩阵逆存在.) B=eye(x); A=A,B; %把A矩阵和单位阵B按行合并 y=2*y; for k=1:x %该for循环是在选主元 max=abs(A(k,k); r=k; for L=k+1:x if max A=unifrnd(0,100,6,6),F=solvenif2(A)%A=unifrnd(0,1

14、00,6,6) 生成0,100区间上的连续型均匀分布6行6列的随机数矩阵A = 95.0129 45.6468 92.1813 41.0270 13.8891 1.5274 23.1139 1.8504 73.8207 89.3650 20.2765 74.6786 60.6843 82.1407 17.6266 5.7891 19.8722 44.5096 48.5982 44.4703 40.5706 35.2868 60.3792 93.1815 89.1299 61.5432 93.5470 81.3166 27.2188 46.5994 76.2097 79.1937 91.6904

15、 0.9861 19.8814 41.8649该矩阵的逆F= F =0.0574 0.0275 0.0365 0.0015 -0.0622 -0.0241 -0.0442 -0.0253 -0.0147 -0.0057 0.0534 0.0156 -0.0139 -0.0061 -0.0211 -0.0009 0.0153 0.0186 -0.0169 -0.0076 -0.0061 -0.0037 0.0313 -0.0060 -0.0364 -0.0461 -0.0471 0.0253 0.0613 0.00900.0272 0.0331 0.0299 -0.0019 -0.0513 -0

16、.0065 A*Fans = 1.0000 -0.0000 0.0000 0.0000 -0.0000 0.0000 0 1.0000 0 -0.0000 0.0000 0.0000 -0.0000 -0.0000 1.0000 -0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0000 1.0000 -0.0000 0 -0.0000 -0.0000 -0.0000 0.0000 1.0000 0.0000 -0.0000 0 0 0.0000 -0.0000 1.00003、利用三角分解的方法（LU分解）编写的求逆算法（代码如下）M函数 函数名：F=solv

17、enif3(A)function L,U,F=solvenif3(A)%实现对矩阵A的LU分解，U为上三角矩阵，L为下三角矩阵，并输出A的逆矩阵n=rank(A);x,y=size(A);if x=y disp(请注意：因为该矩阵非方阵,不能进行LU分解,所以该矩阵无逆.) L=L不存在; U=U不存在; F=A无逆矩阵;else %disp(请注意：该矩阵为方阵.) d=det(A); if d=0 disp(请注意：因为该矩阵非满秩,不能进行LU分解,所以该矩阵无逆.) L=L不存在; U=U不存在; F=A无逆矩阵; else %disp(请注意：因为该矩阵满秩，所以该矩阵逆存在.) B

18、=eye(n); A=A,B; L=zeros(n,n); U=zeros(n,2*n); for i=1:n L(i,i)=1; end for k=1:n for j=k:2*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,k); end end C=U(:,n+1:2*n); U=U(:,1:n); F=UC; end disp(U为上三角矩阵，L为下三角矩阵，F为A的逆矩阵)end注：利用三角分解法进行LU分解具有唯一性。上

19、述的三角分解为直接分解，没有进行选主元的过程，故当A矩阵出现病态时，可能导致结果不正确。对选主元的三角分解这里不在进行详细说明了。有兴趣的读者请看附录一的选主元三角分解求逆矩阵的Matlab代码。运行代码： A=unifrnd(0,100,5,5);L,U,F=solvenif3(A)U为上三角矩阵，L为下三角矩阵，F为A的逆矩阵L = 1.0000 0 0 0 0 0.7267 1.0000 0 0 0 0.8846 1.3486 1.0000 0 0 0.5730 0.9633 1.1734 1.0000 0 0.7428 1.1373 0.1601 -0.6309 1.0000U = 5

20、8.2792 22.5950 20.9069 56.7829 41.5375 0 41.5616 22.7894 38.1588 0.3159 0 0 29.1057 -95.7702 50.2685 0 0 0 103.3674 -81.5894 0 0 0 0 -13.9437F = 0.0307 -0.0596 -0.0781 0.0493 0.0950 -0.0129 -0.0129 -0.0660 0.0284 0.0867 -0.0011 0.0065 0.0532 -0.0075 -0.0624 -0.0043 0.0360 0.0396 -0.0260 -0.0566 -0.0

21、056 0.0381 0.0646 -0.0452 -0.0717 A*Fans = 1.0000 0.0000 0.0000 0 0.0000 -0.0000 1.0000 -0.0000 0.0000 0.0000 0.0000 0.0000 1.0000 0.0000 0.0000 -0.0000 0.0000 0.0000 1.0000 0.00000.0000 -0.0000 -0.0000 0.0000 1.00002.3 三种算法的简单比较运用Matlab7.0中的计时函数对上述三种求逆的算法进行比较：随机生成一个20阶的方阵，分别运用这两种方法求逆，并计算运算时间。 A=uni

22、frnd(0,100,20,20)结果略 tic;F=solvenif1(A);toc该矩阵的逆F= Elapsed time is 0.047000 seconds. tic;F=solvenif2(A);toc该矩阵的逆F= Elapsed time is 0.015000 seconds. tic;L,U,F=solvenif3(A);tocU为上三角矩阵，L为下三角矩阵，F为A的逆矩阵Elapsed time is 0.190000 seconds.仅从计算时间上看，分别运用这三种算法对同一个20阶方阵求逆，初等行变换求逆的效率与另外的两种相比较好。三、特殊矩阵的求逆方法3.1对称正定

23、矩阵求逆设C为n阶对称正定矩阵，y、x为n维向量，其关系式为：y=Cx （3-1）确定了上的一个映像，如能写出逆关系： x=By （3-2）则B为C的逆阵，即B=现将式（3-1）写成 . . （3-3） . .因C对称正定，必有0,用除式（3-3）第一个方程的两端，解出，把和的位置交换，并将代入其他各式得 . . （3-4） . . 事实上式（3-4）可以改写为 . + . + . （3-5） . +如果对式（3-5）中的变量按如下规则重新编号 （3-6）经n次变换后恢复原状，采用变量循环重新编号法的计算公式如下：对于 由于变量循环重新编号法求逆均在下三角阵（包括对角元素）进行，因而运行速快。

24、3.2 有理矩阵求逆定义3.2.1 设其中 对任意的都是有理数，则称A为有理矩阵.引理3.2.1 n阶有理矩阵A与n阶有理矩阵B的乘积可在时间内完成.下面讨论n阶有理矩阵的求逆问题.设，A为有理矩阵，存在.记 其中.若存在，令 由分块矩阵的求逆公式有 显然，都是有理矩阵.令表示n阶有理矩阵求逆的复杂度，应用引理1，计算的运算量不超过，而矩阵加法仅需运算量，于是求的运算量为.用同样的方法讨论计算的运算量可得关系式 令n=，则 显然，从而 若,则有正整数k满足 于是只要将有理矩阵A扩充为阶有理矩阵,其中 且 则求可在时间内完成.由 得 定理3.2.1 n阶有理矩阵的求逆可在时间内完成.四、一般矩阵

25、的广义逆矩阵广义逆矩阵是通常逆矩阵的推广。在许多实际问题中所遇到的矩阵往往是奇异方阵或是任意的矩阵 （一般）, 显然不存在通常的逆矩阵, 这就促使人们去想象能否推广逆的概念。矩阵的几种广义逆对任意复数矩阵, 如果存在复矩阵,满足 (1.1) (1.2) (1.3) (1.4) 则称为的一个 MoorePenrose 广义逆, 并把上面四个方程叫做 MoorePenrose 方程, 简称 MP方程.由于 MP 的四个方程都各有一定的解释, 并且应用起来各有方便之处， 所以出于不同的目的, 常常考虑满足部分方程的 X, 叫做弱逆, 为引用的方便, 我们给出如下的广义逆矩阵的定义.定义4.1 设,

26、若有某个, 满足 MP 方程（1.1）（1.4）中的全部或其中的一部分, 则称为的广义逆矩阵.例如有某个, 只要满足式（1.1） , 则为的广义逆, 记为； 如果另一个, 满足式(1.1), (1.2)则为的广义逆, 记为； 如果, 则同时满足四个方程, 它就是 MoorePenrose广义逆, 等等. 总之, 按照定义 1.1可推得, 满足1个, 2个, 3个, 4个MoorePenrose方程的广义逆矩阵共有种, 但应用较多的是以下五种, , , , .其中每一种广义逆矩阵又都包含着一类矩阵, 分述如下:1: 其中任意一个确定的广义逆, 称作减号逆, 或 逆, 记为； 2: 其中任意一个确

27、定的广义逆, 称作自反广义逆, 记为；3: 其中任意一个确定的广义逆, 称作最小范数广义逆, 记为；4: 其中任意一个确定的广义逆, 称作最小二乘广义逆, 记为；5: 唯一，称作加号逆, 或伪逆, 或 Moore-Penrose 逆, 记为. 为叙述简单起见, 下面我们以 及实矩阵为例进行讨论, 对于及复的矩阵也有相应结果. 本文着重介绍减号逆和加号逆以及左逆与右逆的性质及计算, 并讨论它们在解线性方程组中的应用.4.1 的定义与计算定义4.1.1 设, 若满足, 则称为的记为,由定义可知.例如设, 则就是的, 这里可以任取. 不难看出的逆并不唯一.定理4.1.1 设, , 分别为阶与阶非奇异

28、方阵, 且则 . （证明见参考文献7）例1 求矩阵的广义逆.解 构造分块矩阵, 通过适当变化, 将进行行列变换化为形式, 并求出变换, .,因此有, . 于是我们取, , 均为0得.4.2 加号逆的性质及计算定义4.2.1设, 若存在 阶矩阵 , 它同时满足： 1) 2） 3） 4）则称为 的加号逆, 或伪逆, 或 Moore-Penrose 逆, 记为.从定义中可看出, 加号逆必同时是减号逆、自反广义逆、最小范数广义逆和最小二乘广义逆, 在四个条件中, 与 完全处于对称地位. 因此也是的加号逆, 即有; 另外可见, 加号逆很类似于通常的逆阵, 因为通常的逆也有下列四个类似的性质: 1. 2.

29、 3. 4. 由定义1.2.1 中的条件 3）和 4）还可看出, 与都是对称矩阵.前面已经介绍了什么样的矩阵称为广义逆矩阵, 下面将讨论广义逆矩阵的唯一性.定理4.2.1 对任意, 存在且唯一.证明 设, 若则是阶零矩阵, 显然阶零矩阵满足条件. 若则的满秩分解为, 其中, , 于是即为所求的.因为 (1) ；(2) ；(3) ；(4) .由此说明了广义逆的存在性.又设则有 .这便说明了的唯一性.定理4.2.2 设为秩为的矩阵, 其满秩分解为, 其中, , 则.的唯一性前面已经作出了说明, 此定理的证明见参考文献74.3 左逆与右逆的定义定义4.3.1 设是矩阵, 若有矩阵满足(或), 则称为

30、的右逆(或左逆), 记为(或).定理4.3.1 设是的矩阵, 有右(左)逆()的充要条件是().若有右(左)逆, 则其中一个右(左)逆是(), 通式为()其中是任意满足的矩阵.证明 充分性: 已知, 则, 是可逆矩阵, 若记, 则, 因此是的右逆.必要性: 设是的一个右逆, 则. 由于,因此.设是任意满足的矩阵, 最后证明右逆的通式可以表示成为的形式.由于, 因此是的右逆. 设是的任意右逆,记, 则因此. 又因为=,由上分析可知的任意右逆都可找到使其表示为的形式.因此矩阵的右逆的通式为.对于左逆同理证明.例 求矩阵的左逆.解 由于,所以我们有例 设 ,试求其右逆.解 易知rank,即是最大秩矩

31、阵，有 =.结束语矩阵理论不但是经典数学的基础，同时又是很有实用价值的数学理论。计算机的广泛应用为矩阵理论的应用开辟了广阔的应用前景。逆矩阵的概念在矩阵理论中占有重要位置。随着科技的发展，逆矩阵的应用范围也越来越广,可见对逆矩阵的研究具有深远的意义。本文在逆矩阵的性质以及矩阵可逆的等价条件的基础上介绍了定义法、伴随矩阵法、初等变换法、三角分解法、分块矩阵法等求逆矩阵的常用方法。但由于矩阵的类型繁多，对于不同类型的矩阵可以用不同的方法来求，这样可以达到简单、易求的目的。附录一 M函数名 L,U,F=solvenif4(A)function L,U,F=solvenif4(A)%实现对矩阵A的LU

32、分解，U为上三角矩阵，L为下三角矩阵，并输出A的逆矩阵n=rank(A);x,y=size(A);if x=y disp(请注意：因为该矩阵非方阵,不能进行LU分解,所以该矩阵无逆.) L=L不存在; U=U不存在; F=A无逆矩阵;else %disp(请注意：该矩阵为方阵.) d=det(A); if d=0 disp(请注意：因为该矩阵非满秩,不能进行LU分解,所以该矩阵无逆.) L=L不存在; U=U不存在; F=A无逆矩阵; else %disp(请注意：因为该矩阵满秩，所以该矩阵逆存在.) B=eye(n); A=A,B; L=zeros(n,n); U=zeros(n,2*n); p=eye(n);%p记录了选择主元时候所进行的行变换 for k=1:n-1 r,m=max(abs(A(k:n,k); %选列主元 m=m+k-1; if(A(m,k)=0) if(m=k) A(k