大型稀疏矩阵的LU分解及特征值求解

稀疏矩阵方程算法稀疏矩阵是指矩阵中绝大多数元素为0的矩阵。

在实际问题中，很多矩阵都是稀疏的，例如图像处理、自然语言处理等领域。

由于稀疏矩阵的特殊性，传统的矩阵运算方法效率较低，因此需要设计高效的算法来解决稀疏矩阵方程。

稀疏矩阵方程是指形如Ax=b的线性方程，其中A是一个稀疏矩阵，b是一个向量。

解决稀疏矩阵方程的一种常用方法是使用迭代算法，例如共轭梯度法（Conjugate Gradient，CG）和广义最小残差法（Generalized Minimal Residual，GMRES）等。

共轭梯度法是一种迭代法，它可以用来解决对称正定稀疏矩阵方程。

该方法的基本思想是通过最小化残差的二次范数来逼近方程的解。

具体而言，共轭梯度法通过迭代计算一个与残差正交的搜索方向，并在该方向上进行搜索，直到找到方程的解。

广义最小残差法是一种迭代法，它可以用来解决非对称稀疏矩阵方程。

该方法的基本思想是通过最小化残差的二范数来逼近方程的解。

与共轭梯度法不同的是，广义最小残差法使用Krylov子空间来进行搜索，并在该子空间上进行最小化残差的计算。

除了迭代算法外，还可以使用直接求解方法来解决稀疏矩阵方程。

其中一种常用的方法是LU分解。

LU分解是将稀疏矩阵A分解为一个下三角矩阵L和一个上三角矩阵U的乘积，即A=LU。

通过LU分解，可以将原始方程Ax=b转化为Ly=b和Ux=y两个方程，进而求解出x的值。

稀疏矩阵方程的求解算法还有很多，例如Jacobi迭代法、高斯-赛德尔迭代法等。

这些算法在不同的问题和应用中具有不同的优势和适用性。

在实际应用中，稀疏矩阵方程的求解是一个复杂且关键的问题。

通过选择合适的算法和优化技术，可以提高求解的效率和精度。

同时，还可以利用稀疏矩阵的特殊性质，例如压缩存储和并行计算等，进一步提高算法的性能。

稀疏矩阵方程是一类特殊的线性方程，传统的矩阵运算方法在处理稀疏矩阵时效率较低。

针对稀疏矩阵方程，可以采用迭代算法和直接求解方法来求解。

矩阵的LU分解是一种将一个矩阵分解为一个下三角矩阵（L）和一个上三角矩阵（U）的乘积的方法。

LU分解是线性代数中一种重要的矩阵分解方法，它在求解线性方程组、计算行列式、计算矩阵的逆等方面都有广泛的应用。

对于一个$n \times n$的矩阵A，其LU分解的计算量主要取决于以下几个因素：
1. 存储量：LU分解需要存储原始矩阵A、下三角矩阵L和上三角矩阵U，因此需要额外的$3n^2$个浮点数存储空间。

2. 计算量：LU分解需要进行一系列的矩阵乘法和行交换操作，因此计算量相对较大。

具体来说，LU分解的计算量主要包括以下几部分：
* 计算L矩阵：需要执行$n(n+1)/2$次乘法操作，其中$n$是矩阵的阶数。

* 计算U矩阵：需要执行$n^3/6$次乘法操作，其中$n$是矩阵的阶数。

* 行交换操作：需要执行$n-1$次行交换操作，其中$n$是矩阵的阶数。

因此，LU分解的计算量大约为$O(n^3)$，其中$n$是矩阵的阶数。

在实际应用中，为了提高计算效率，通常会采用一些优化算法和并行计算技术来加速LU分解的计算过程。

lu分解方法求矩阵特征值
矩阵特征值可以通过LU分解方法来求解。

首先，我们需要将矩阵A进行LU分解，得到一个下三角矩阵L和一个上三角矩阵U，使得LU=A。

接下来，我们可以利用LU分解的结果来求解特征值。

假设我们已经得到了LU分解后的矩阵L和U，那么我们可以进行如下步骤来求解矩阵A的特征值：
1. 首先，求解U矩阵的特征值。

U矩阵是一个上三角矩阵，其特征值就是它的对角线元素。

2. 接下来，求解L矩阵的特征值。

L矩阵是一个下三角矩阵，其特征值也是它的对角线元素。

3. 最后，将U矩阵和L矩阵的特征值组合起来，就得到了矩阵A的特征值。

需要注意的是，LU分解方法求解特征值的过程相对比较复杂，尤其是对于大型矩阵而言。

在实际应用中，可以借助计算机软件来进行LU分解和特征值求解，以提高计算的准确性和效率。

总之，通过LU分解方法可以求解矩阵的特征值，但需要注意计算的复杂性和精度的要求。

希望这个回答能够帮助到你理解如何利用LU分解方法来求解矩阵的特征值。

lu分解算法
LU分解算法是一种将一个非奇异矩阵分解成一个下三角矩阵L和一个上三角矩阵U的方法，它可以用于解线性方程组以及求矩阵的逆等计算中。

具体的LU分解算法如下：
输入：一个n×n的非奇异矩阵A
输出：下三角矩阵L和上三角矩阵U
1. 初始化一个n×n的下三角矩阵L和一个n×n的上三角矩阵U，使它们的所有对角元素为1。

2. 对于矩阵A的第一行，将其作为矩阵U的第一行。

3. 对于矩阵A的第一列，将其除以矩阵U的第一个元素得到矩阵L的第一列。

4. 对于矩阵A的剩余行，以及对应的列，进行如下操作：
- 计算当前元素的值，即A(i, j)减去矩阵L的第i行与矩阵U的第j列的内积。

- 如果i小于等于j，将计算得到的值赋给矩阵U的第i行第j列元素。

- 如果i大于j，将计算得到的值除以矩阵U的第j列第j个元素，然后赋给矩阵L的第i行第j列元素。

5. 返回矩阵L和矩阵U作为结果。

通过LU分解算法，可以将解线性方程组的计算转化为简单的矩阵乘法和求解步骤。

此外，通过求解LU分解后的矩阵，还可以求矩阵的逆和行列式等相关计算。

大型稀疏矩阵直接求解算法的研究及实现共3篇大型稀疏矩阵直接求解算法的研究及实现1大型稀疏矩阵直接求解算法的研究及实现随着计算机技术的不断发展和数学建模需求的增加，大型稀疏矩阵直接求解算法的研究和实现日益受到人们的关注。

在实际应用中，大型稀疏矩阵经常出现在各种科学计算、工程计算以及机器学习等领域。

因此，如何高效地求解大型稀疏矩阵成为了一个十分重要的问题。

一般来说，大型稠密矩阵的求解可以使用各种经典算法，如高斯消元、LU分解等。

然而，大型稀疏矩阵的求解却需要特殊的算法和数据结构。

传统的直接求解方法存在着效率低下和存储空间过大等问题，因此研究者们提出了许多改进方法和优化方案。

稀疏矩阵存储结构是求解算法中的重要问题之一。

目前，广泛应用的稀疏矩阵存储格式包括压缩列（Compressed Column，CC）、压缩行（Compressed Row，CR）以及双重压缩（Double Compressed）等。

这些存储格式各有优缺点，具体用哪一种存储格式取决于矩阵的具体特点和求解算法的需求。

比如，在随机梯度下降等机器学习算法中，常常使用压缩行存储方式来优化矩阵乘法操作的速度。

多核并行、GPU加速等技术也被广泛应用于大型稀疏矩阵的求解算法中，以提高计算效率。

并行求解算法可以将巨大的计算任务划分成多个子任务，并分配给多个核心同时执行，充分利用计算机的计算资源。

而GPU加速则充分利用了GPU的特殊架构，通过将计算任务映射到各个流处理器上并行执行，进一步提高求解效率。

除了以上所述的算法优化和技术应用，近年来还出现了一些新的求解算法。

比如，基于埃米尔特矩阵分解的求解算法，具有比传统LU分解更快的求解速度；基于内点法的求解算法，在高稀疏性的情况下，具有比其他算法更优的求解速度和精度。

综上所述，大型稀疏矩阵直接求解算法的研究和实现是一个充满挑战的领域。

在实际应用中，选择适合的算法和存储结构，并结合多核并行、GPU加速等技术，可以有效提高求解速度和精度。

矩阵LU分解的几种算法Doolittle分解将矩阵A分解为单位下三角矩阵L和上三角矩阵UCrout 分解将矩阵A分解为下三角矩阵L和单位上三角矩阵UCholesky分解Doolittle分解和Crout 分解适于一般非奇异的矩阵，但对于一些更特殊的矩阵，我们有更好的分解方法。

基础概念矩阵A对称：A^T=A矩阵A正定：A的各阶顺序主子式大于0，对于实对称矩阵A正定的等价条件是A的特征值全为正假设矩阵A是对称正定矩阵，则可以分解为：其中L为下三角矩阵注：这里不给出证明，具体的分解过程，大部分数学软件都有相应的函数，我们更关心如何应用这样可以将求解线性方程组的过程看做两个步骤由于L为下三角矩阵，所以x,y都很好求解，简化了运算过程。

现在假设A为对称矩阵，去掉正定的条件，但是规定矩阵A的各阶顺序主子式不为0那么矩阵A可以做如下分解其中D为对角阵，L为下三角矩阵这样我们可以将求解线性方程组的过程同样看做两个步骤由于D为对角阵，它的逆就是它的倒数，其余的矩阵都是三角矩阵，所以计算也十分简便。

值得注意的是，显然，如果矩阵A是对称正定的，那么也是可以分解为LDL^T的，但如果矩阵A 不是正定的，那么不能分解为LL^T。

补充知识：一个三角矩阵的逆，也是三角矩阵且对角线上元素是倒数关系，但其余位置不是的。

例如：追赶法追赶法是针对带状矩阵（尤其是三对角矩阵）这一大稀疏矩阵的特殊结构，得出的一种保带性分解的公式推导，实质结果也是LU分解可以将一个三对角的稀疏矩阵分解为如下形式：其中三对角矩阵A为：最后提及一句：mathematica中提供LUDecomposition,CholeskyDecomposition 两个函数实现矩阵的LU分解。

实验二矩阵的LU分解一、题目：求矩阵[2 1 1 2;1 2 3 2;2 4 1 1;3 1 2 3]的Doolittle分解二、方法：Doolittle分解，矩阵的紧凑格式的LU分解法三、程序：(1) 紧凑格式的LU分解法function a=nalupad(a)n=length(a);a(2:n-1)=a(2:n-1)/a(1,1);for k=2:na(k,k:n)=a(k,k:n)-a(k,1:k-1)*a(1:k-1,k:n);a(k+1:n,k)=(a(k+1:n,k)-a(k+1:n,1:k-1)*a(1:k-1,k))/a(k,k);end(2) LU分解function f=LU_decom(A)[m,n]=size(A)L=eye(n);U=zeros(n);flag='ok';for i=1:nU(1,i)=A(1,i);endfor r=2:nL(r,1)=A(r,1)/U(1,1);endfor i=2:nfor j=i:nz=0;for r=1:i-1z=z+L(i,r)*U(r,j);endU(i,j)=A(i,j)-z;endif abs(U(i,i))<epsflag='failure'return;endfor k=i+1:nm=0;for q=1:i-1m=m+L(k,q)*U(q,i);endL(k,i)=(A(k,i)-m)/U(i,i);endendLUend四、结果:(1) 紧凑格式的LU分解法>> nalupad([2 1 1 2;1 2 3 2;2 4 1 1;3 1 2 3])ans =2.0000 1.0000 1.0000 2.00000.5000 1.5000 2.5000 1.00001.00002.0000 -5.0000 -3.00003.0000 -1.3333 -0.4667 -3.0667(2) LU分解>> LU_decom([2 1 1 2;1 2 3 2;2 4 1 1;3 1 2 3])m =4n =4L =1.0000 0 0 00.5000 1.0000 0 01.00002.0000 1.0000 01.5000 -0.3333 -0.2667 1.0000U =2.0000 1.0000 1.0000 2.00000 1.5000 2.5000 1.00000 0 -5.0000 -3.00000 0 0 -0.4667五、拓展：矩阵分解成LU形式是有条件的，首先矩阵必须是非奇异的矩阵，其次矩阵的全部顺序主子式非零的时候才能完全保证矩阵可分解成LU且分解唯一。

浅谈矩阵的LU 分解和QR 分解及其应用基于理论研究和计算的需要，往往有必要把矩阵分解为具有某种特性的矩阵之积，这就是我们所说的矩阵分解.本文将介绍两种常用的矩阵分解方法，以及其在解线性方程组及求矩阵特征值中的应用.1.矩阵的LU 分解及其在解线性方程组中的应用 1.1 高斯消元法通过学习，我们了解到利用Gauss 消去法及其一些变形是解决低阶稠密矩阵方程组的有效方法.并且近些年来利用此类方法求具有较大型稀疏矩阵也取得了较大进展.下面我们就通过介绍Gauss 消去法，从而引出矩阵的LU 分解及讨论其对解线性方程组的优越性. 首先通过一个例子引入：例1,解方程组(1.1)(1. 2)(1.3)解. 1Step (1.1)(2)(1.3)⨯-+ 消去(1.3)中未知数,得到23411x x --=- (1.4)2Shep . (1.2)(1.4)+ 消去(1.4)中的未知数2x有12323364526x x x x x x ++=-=-=-⎧⎪⎨⎪⎩ 显然方程组的解为*x =123⎛⎫ ⎪ ⎪ ⎪⎝⎭上述过程相当于 111604152211⎛⎫ ⎪- ⎪ ⎪-⎝⎭~111604150411⎛⎫ ⎪- ⎪ ⎪---⎝⎭~111604150026⎛⎫⎪- ⎪ ⎪--⎝⎭2-（）+ ()i i r 表示矩阵的行由此看出，消去法的基本思想是：用逐次消去未知数的方法把原方程化为与其等价的三角方程组.下面介绍解一般n 阶线性方程组的Gauss 消去法.设111n n1nn a a a a A ⎛⎫ ⎪=⎪ ⎪⎝⎭ 1n x X x ⎛⎫ ⎪= ⎪ ⎪⎝⎭ 1n b b b ⎛⎫ ⎪= ⎪ ⎪⎝⎭则n 阶线性方程组AX b = (1.5)并且A 为非奇异矩阵.通过归纳法可以将AX b =化为与其等价的三角形方程，事实上： 及方程(1.5)为()()11A X b =,其中()1A A = ()1b b =(1) 设(1)110a≠，首先对行计算乘数()()11i1111i a m m =.用1i m -乘(1.5)的第一个方程加到第()2,3,,i i n =⋯个方程上.消去方程(1.5)的第2个方程直到第n 个方程的未知数1x .得到与(1.5)等价的方程组()()()11n 12n 111nn 0a a x x a ⎛⎫⋯⎛⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⋯⎝⎭⎝⎭=()()112n b b ⎛⎫ ⎪⎪ ⎪⎝⎭简记作()()22Ab = (1.6)其中()()()()()()211211111 ijij i ij i i i a m b b m a a b =-=- (2) 一般第()11k k n ≤≤-次消去，设第1k -步计算完成.即等价于()()k k AX b = (1.7)且消去未知数121,,,k x x x -⋯.其中()()()()()()()()()()1111112122222k k k k kk knk nknna n n a a a a a A a a a a ⎛⎫⎪ ⎪⎪ ⎪= ⎪ ⎪ ⎪ ⎪⎝⎭设()0k kk a ≠计算()() (i=/1,,)k k ik ikkkaa k m n =+⋯,用()()(1,,)a n ikik k kki k n a m ==+⋯消去第1k +个方程直到第n 个方程的未知数k x .得到与(1.7)等价的方程组()()1k 1k A X b ++= 故由数学归纳法知，最后可以把原方程化成一个与原方程等价的三角方程组.但是以上分析明显存在一个问题，即使A 非奇异也无法保证()0i ii a ≠,需要把非奇异的条件加强.引理1 约化主元素()01,,)i ii a k ≠=⋯（i 的充要条件是矩阵A 的顺序主子式0i D ≠.即1111110,0ikk kkk a a D a a D a =≠=≠⋯证明 利用数学归纳法证明引理的充分性.显然，当1k = 时引理的充分性是成立的，现在假设引理对1k -是成立的，求证引理对k 亦成立.有归纳法，设()()01,21iii i a k ≠=⋯-于是可用Gauss 消去法将中，即()()()()()()()()()()()11111121n22222n 1k k k k k kk knnknn a a a a a A a a a a A ⎛⎫ ⎪ ⎪⎪ ⎪→= ⎪ ⎪ ⎪ ⎪⎝⎭即()()()()()()()()11121231112211223112233222a a D a D a a a a a ===()()()()()()()()()11111122212222k 11122k k k kkk kka a a a a a D a a a ==⋯ (1.8) 由设0(1,,)i i D k ≠=⋯及式(1.8)有()0k kk a ≠显然，由假设()()01,2iiii k a ≠=⋯，利用(1.8)亦可以推出0(1,,)i i D k ≠=⋯ 从而由此前的分析易得；定理1 如果n 阶矩阵A 的所有顺序主子式均不为零，则可通过Gauss 消去法（不进行交换两行的初等变换），将方程组(1.5)约化成上三角方程组，即()()()()()()()()()1111111121122222222b b b n n n n nn n n a a a x x a a x a ⎛⎫⎛⎫⎛⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪= ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎪⎝⎭⎝⎭⎝⎭ (1.9) 1.2 矩阵LU 分解从而由以上讨论即能引出矩阵的LU 分解，通过高等代数我们得知对A 施行行初等变换相当于用初等矩阵左乘A ，即()()()()121211L A Lb A b == 其中 211n11101L m m ⎛⎫⎪- ⎪= ⎪⎪-⎝⎭一般第k 步消元，，相当于()()()()11k kk k k kL A A L b b ++==重复这一过程，最后得到()()()()11211121n n n n Ab L L L A L L L b --⎧⋯=⎪⎨⋯=⎪⎩ (1.10) 其中k 1,111m 1n k k km L +⎛⎫ ⎪ ⎪ ⎪=⎪- ⎪ ⎪ ⎪ ⎪-⎝⎭将上三角形矩阵()n A U 记作,由式(1.9)得到111121=U n A L L L LU ----⋯=,其中211111211211m 1n n n m L L m L L ----⎛⎫⎪⎪=⋯= ⎪⎪⎝⎭由以上分析得；定理2 （LU 分解） 设A 为n 阶矩阵，如果A 的顺序主子式i 0(1,2,,1)D i n ≠=-.则A 可分解为一个单位下三角矩阵L 和一个上三角矩阵U的乘积，且这种分解是唯一的.证明 由先前的分析得出存在性是显然的，即A LU =.下证唯一性,设A LU CD == 其中L , C 为单位下三角矩阵，U ,D 为上三角矩阵.由于1D -11D C L U --=上式右端为上三角矩阵，左端为单位下三角矩阵，从而上式两端都必须等于单位矩阵，故U D =,L C =.证毕.例2 对于例子1 系数矩阵矩阵111041221A ⎛⎫ ⎪=- ⎪ ⎪-⎝⎭由Gauss 消去法，得结合例1，故100111010041211002A LU ⎛⎫⎛⎫⎪⎪==- ⎪⎪ ⎪⎪--⎝⎭⎝⎭对于一般的非奇异矩阵，我们可以利用初等排列矩阵kki I (由交换单位矩阵I的第k 行与第k i 行得到)，即()()()()()()()()111212111111,,kk k k k ki k i k k i i k A L b L I A I b L I A I b A L b++⎧==⎪⎨==⎪⎩ (1.11) 利用(1.11)得()1111,11n nn n i i L I L U I A A ---==.简记做.其中下面就n 情况来考察一下矩阵.()()4321444343544332211443443243)(i i i i i i i i i i I L I L I L I A I L I I L I I A A L L I ===⨯4324324321432()i i i i i i I I I L I I I 43214321 )(i i i i I I I I A从而记从而容易的为单位下三角矩阵,总结以上讨论可得如下定理.定理3 如果A 非奇异矩阵，则存在排列矩阵P 使PA LU = 其中L 为单位下三角矩阵，U 为上三角矩阵.1.3 矩阵LU 分解的应用以上对非奇异矩阵A 的LU 分解进行了全面的讨论，一下我们就简单介绍一下应用.对于矩阵A 一旦实现了LU 分解，则解线性方程的问题,便可以等价于：(1)Ly b = 求y (2)=Ux y , 求x (1.12)即，设A 为非奇异矩阵，且有分解式A LU =，其中L 为单位下三角矩阵，U 为上三角矩阵。

增量法求解大规模稀疏矩阵方程组近年来，随着人工智能、大数据分析和科学计算等领域的快速发展，大规模稀疏矩阵方程组的求解变得尤为重要。

在实际应用中，由于矩阵规模大、非零元素分布稀疏等特点，传统的直接或间接法求解效率低下，甚至无法胜任。

增量法成为了求解大规模稀疏矩阵方程组的一种有效手段。

在本篇文章中，我们将从简入深地介绍增量法的基本原理、常见算法以及应用实例，帮助读者全面、深刻地理解增量法在求解大规模稀疏矩阵方程组中的重要作用。

1. 基本原理大规模稀疏矩阵方程组中的“大规模”指的是矩阵的规模非常大，通常是上万甚至上亿级别；“稀疏”则表示矩阵的非零元素分布非常稀疏，大部分元素为零。

在实际问题中，这样的矩阵可能来自于网络分析、有限元法、信号处理等领域。

传统的直接法求解（如高斯消元法）由于需要处理大量的零元素而效率较低，而间接法求解（如迭代法）则常常陷入收敛速度慢、计算精度不高等问题。

增量法的基本原理是将原始的大规模稀疏矩阵方程组逐步化为一个个规模较小的子问题，并利用迭代或直接法求解这些子问题，最终得到原问题的解。

通过这种方式，增量法在一定程度上克服了传统方法的缺点，提高了求解效率和精度。

2. 常见算法在增量法求解大规模稀疏矩阵方程组中，有几种常见的算法被广泛应用，包括但不限于以下几种：2.1 LU分解LU分解是一种经典的增量法求解大规模稀疏矩阵方程组的方法。

它将原始矩阵分解为一个下三角矩阵和一个上三角矩阵的乘积形式，从而将原问题转化为两个规模较小的子问题。

可以通过前代和回代的方式分别求解这两个子问题，从而得到原问题的解。

2.2 Cholesky分解Cholesky分解是一种针对对称正定矩阵的增量法求解方法。

与LU分解类似，Cholesky分解也是将原始矩阵分解为一个下三角矩阵和其转置矩阵的乘积形式。

通过这种分解，可以将原问题转化为一个规模较小的子问题，并利用直接法求解得到原问题的解。

2.3 前代和回代前代和回代是一种常见的增量法求解大规模稀疏矩阵方程组的迭代方法。

几类稀疏矩阵特征值的算法实现吴承逊;谌稳【摘要】总结了几类适用于用迭代法求解稀疏矩阵特征值的算法.文章用到了幂法和反幂法,并在幂法算法的基础上对其进行规范化,说明了二分法求特征值的原理;实现了用二分法求对称三对角矩阵所有特征值的算法;讨论了能把大部分稀疏矩阵变成对称三对角矩阵的Lanczos方法,且对不同类稀疏矩阵使用Lanczos方法进行探讨,把二分法和Lanczos方法结合到一个算法中;并通过数值实验验证了这些算法的有效性.【期刊名称】《湖南城市学院学报（自然科学版）》【年(卷),期】2018(027)005【总页数】5页(P51-55)【关键词】稀疏矩阵;特征值;幂法;Lanczos法;算法实现【作者】吴承逊;谌稳【作者单位】湖南城市学院理学院,湖南益阳 413000;湖南城市学院理学院,湖南益阳 413000【正文语种】中文【中图分类】O242迭代法可以用来近似解决稀疏矩阵的特征值问题，它有存储空间小、程序简单的特点，是求解稀疏矩阵特征值的有效方法，其主要方法包括幂法及反幂法、二分法和Lanczos方法．幂法及反幂法是求稀疏矩阵特定特征值的常用方法，它适用性强，但是不能求出全部的特征值．随着J. W. Givens和A. S. Householder所研究的Givens变换和Householder变换的出现，使得任何实对称矩阵都可以通过这2种正交变换约化变成对称的三对角矩阵，于是二分法[1]得以出现，它的理论基础是Sturm序列．二分法是解决对称的三对角矩阵全部特征向量和特征值的有效方法．但是由于稀疏矩阵的阶数过于庞大，约化矩阵变得十分困难和复杂，且二分法不能很好的求解非对称矩阵．Lanczos[2]研究出了一种Lanczos方法，把稀疏矩阵化成对称三对角矩阵，这样稀疏矩阵的特征值问题就变成了转化的对称三对角矩阵的特征值问题，然后用QR方法或者二分法很容易就求解出来了．然而，实际计算中Lanczos向量很容易失去正交性，且长期是不稳定的．1970年之后，Peters Wilkinson、Golub等人重新研究了该方法，对Lanczos方法进行改进，于是该方法成为能够求解大型的稀疏矩阵特征值问题的有效方法[3]．近年来，关于Lanczos方法的研究主要在其应用方面，曾玮、赵永华[4]基于谱分割并行求解稀疏矩阵特征值，并将矩阵的特征值求解区间划分为多个独立的子区间，分别对各个子区间内的特征值进行独立求解﹒焦宝宝[5]用重正交化Lanczos法求解大型非正交归问题，计算结果表明了Lanczos法在Matlab编程和大型壳模型计算中的精确性和可行性﹒幂法是用来求解矩阵的主特征值及其对应的特征向量的迭代方法，其基本思路是：设Ａ是一个有着n个线性无关特征向量的实矩阵，于是它的特征值是，对应的特征向量是．假设Ａ的主特征值是实根，且有，任取1个非零初始向量，用A构造的一个向量序列为可称作是迭代向量．幂法算法的步骤如下：(1)初始化设置，随机生成稀疏矩阵A，是第1个元素为1，其余为0的且和A矩阵同阶的列向量；(2)计算，；(3)计算的绝对值，如果满足误差精度，就输出主特征值m；否则，下标加1，回到步骤(2)．运用幂法求解的主特征值m和它对应的特征向量v，若m的绝对值大于1或者小于1，就会使得迭代向量v的每个不等于０的分量随着k趋近于无穷大时而趋向于无穷大(或者趋向于0)，这样在电脑上用算法去实现时就可能会“溢出”，为了解决这个缺点，需要将幂法中的迭代向量进行规范化．对于迭代向量，将它规范化，可得到其中，max{}是指v的绝对值最大的分量．于是，规范化幂法就是按照公式(3)进行构造的向量序列，即由式(3)得到的最大特征值的近似值为，它对应的近似特征向量的迭代向量为．规范化幂法算法的步骤如下：(1)初始化设置，随机生成稀疏矩阵A，是第1个元素为1，其余为0的且和A矩阵同阶的列向量；(2)计算，；(3)计算的绝对值，如果满足误差精度，就输出主特征值；否则，下标加1，回到步骤(2)．反幂法是幂法的一种变形，它用来求解矩阵模最小的特征值及其对应的特征向量，也能用于计算给定的近似特征值的特征向量．幂法是求A的最大特征值的问题，而反幂法则是先求解的最大特征值，进而求得A的最小特征值的问题．反幂法的迭代公式可描述为：对于任意的初始向量，可以构造向量序列，即其中，迭代向量能够求解线性方程组得出；最小特征值为．当然，反幂法是可以用于计算出给定特征值的特征向量的．例如给定A的1个特征值a，其对应的特征向量为x，则可以对(p约等于a)应用反幂法，它的1个特征值是，这个值是所有特征值中最小的，于是就能够求出这个特征向量．同样，由于特征值值小，迭代次数往往也不大，几次就能迭代出特征向量．反幂法求最小特征值与幂法算法相似，只是在迭代公式上有些许不同，用反幂法求给定特征近似值的特征向量步骤如下：(1)初始化设置，输入稀疏矩阵A和给定的特征近似值p，是第1个元素为1，其余为0的且和A矩阵同阶的列向量；(2)计算出，以及；(3)求出，tol是误差精度，如果tol满足规定的条件，则上述循环结束，输出所求的；否则，下标加1，回到步骤(2)．此算法设计并非最优，因为算法中涉及到求矩阵的逆，可以对其进行如下改进：先把A进行LU分解，再通过求解，来进行算法的迭代．幂法及反幂法只能够求出可对角化矩阵A的最大和最小特征值，无法求出A的所有特征值，于是需要用到二分法．二分法的理论基础是Sturm序列，下面介绍下什么是Sturm序列．对于n阶对称的三对角矩阵用来表示的k阶顺序主子式，即令，且满足式(5)的关系：则称是B的Sturm序列，同时从公式(5)中很明显的看出来，是B的特征多项式．把中每相邻2元素符号相同的次数称为同号数，记为．那么，应用二分法的步骤就可以分为以下3步：(1)用Householder变换和Givens变换把所求矩阵约化为对称三对角矩阵；(2)再求出三对角矩阵的顺序主子式符号，用来确定中根的个数，式中分别表示A 的r阶顺序主子阵和单位矩阵的r阶主子，然后再用区间的二分法隔离出所需要的根；(3)计算所求矩阵的相应特征向量，可以通过反幂法求出．二分法的关键是可以根据Sturm序列的同号数来判断特定区间的根的个数，也就是说就是在区间中根的个数，即B中不小于的特征值的个数．特征值问题就是求解式(6)的问题，即其中，A是一个非对称实矩阵；、分别是A的特征值所对应的右特征向量和左特征向量．Lanczos方法就是生成了2组双正交向量和，使得，其中是指n阶的单位矩阵，满足，且由上面2个公式可以得出将式(7)以列的形式描述出来，得到对其进行重新排序，得到递推公式当然，系数不是唯一的，只要两者之一满足即可，再引入双正交的参数，取、以及剩下的系数，其中可以由方程得到．上述Lanczos方法有着显著的优势，它把非对称实矩阵约化成了对称三对角矩阵，上面描述的迭代过程最多进行n次就终止了，除了元素间可能出现的正负号之差外．本文对Lanczos方法的算法[6]，作出一点改进，把它和二分法结合在一起，设计出一个新算法，它们两者是承上启下的，这个算法的优势在于最后输出结果是对称三对角阵、最大和最小特征值，具体算法步骤如下：(1)初始化设置，，；(2)从循环计算，且循环次数为n次(n为矩阵的阶数)．这样，矩阵就被上述循环算法化成了对称的三对角矩阵；(3)对转化的三对角矩阵使用二分法，然后按照公式(5)计算和同号数；再判断是否，如果是，则，否则；继续判断的值是否不小于误差，若是，则输出，否则就重新进行这个循环．可以使用Matlab随机生成稀疏矩阵，如生成1个100阶的稀疏矩阵，为了保证这个稀疏矩阵是一个非奇异矩阵，随机生成100阶矩阵后，取这个矩阵加上1个100阶的单位矩阵的和矩阵．随机选取4组不同阶数的稀疏矩阵来进行幂法的算法运算，稀疏矩阵的非0元素分布密度为0.3，运算结果如表1所示．表1的二范数q是的二范数，其中，A是所求的稀疏矩阵；是A的特征值；是对应的特征向量．显而易见，规范化幂法的算法是更优于幂法算法的，其改进程度也十分显著．从表1中可看出，幂法的缺点十分明显，它在前3组的稀疏矩阵迭代中迭代次数比规范化幂法的迭代次数多，差向量二范数也十分的大，说明它的特征向量中元素的值很大，这样计算会出现舍入误差，造成数值的不准确；在第4组稀疏矩阵的迭代中，它的特征值是1，明显是错误的，因为特征值为1，说明在特别稀疏的情况下，如果主特征值是很小的一个数，趋近于零，那么幂法就会失败．然而，规范化幂法就没有这样的缺点，它是一个更优于幂法的改进算法，能够适用于绝大多数的具有完备特征向量组的稀疏矩阵．随机选出了3个的10阶稀疏矩阵，用Lanczos算法来运算，运算结果如表2所示．从表2的运算情况来看，这个Lanczos算法较为完美的把稀疏矩阵化成了对称三对角矩阵，再通过把最大及最小值用幂法和反幂法验证，发现相差不大，证明了该方法有效．稀疏矩阵的特征值问题长期是矩阵研究中的一项重要任务，求解稀疏矩阵特征值的方法有很多是不稳定的、有误差的，虽然有着各种近似迭代方法的出现和发展，但是在每一种方法都有对应的特定情况，所以要证明哪种方法是最行之有效的，这是十分困难的．本文研究了幂法、规范化幂法、反幂法、二分法和Lanczos方法，发现能够解决稀疏矩阵特征值问题的主要方法还是Lanczos方法；而且，当用Lanczos方法运算后，用二分法求出的最大及最小特征值，有时候会和幂法与反幂法求出的特征值有一些误差．这是因为Lanczos方法迭代的双正交向量在多次迭代后会出现失去正交性的问题，求出的特征值不太理想，这个问题还有待进一步的研究和探讨．【相关文献】[1]黄云清, 舒适, 陈艳萍, 等. 数值计算方法[M]. 北京: 科学出版社, 2009.[2]LANCZOS C. An iteration method for the solution of the eigenvalue problem of linear differential and intergral operators[J]. Journal of Research of National Bureau of Standards, 1950, 45(4): 255-282.[3]周树荃. 计算大型稀疏矩阵部分特征值问题的有效方法[J]. 南京航空学院学报, 1981(1): 1-18.[4]曾玮, 赵永华. 基于谱分割的稀疏矩阵特征值问题并行求解[J]. 数值计算与计算机应用, 2015,36(2): 132-146.[5]焦宝宝. 用重正交化Lanczos法求解大型非正交归一基稀疏矩阵的特征值问题[J]. 物理学报, 2016, 65(19): (192101)1-8.[6]邓健新. 用Lanczos方法解高阶稀疏矩阵广义特征值问题的某些经验[J]. 数值计算与计算机应用, 1980, 1(3): 173-180.。