python 逆矩阵_如何用Python快速求解正定对称矩阵的逆？

正定对称矩阵是一类比较特殊的矩阵。其正定性决定了它的特征值全为正，从而它必然是非奇异的，也就是一定有逆矩阵存在。其对称性使得它可以进行对称分解，从而在进行各种操作时可以有各种便捷的方法选用。

这里我们主要探讨一下对于一个严格的对称正定矩阵，在Python的库里面如何快速求解。

这里我们主要讨论scipy库中的相关方法。scipy是python中矩阵操作应用最为广泛的库之一，它的底层实际上是使用的大名鼎鼎的LAPACK库，其运行效率和计算精度都十分高。众所周知，最专业的数学软件Matlab其实底层也是直接调用的LAPACK库。因此直接从scipy入手，其正确性是十分有保障的。

废话不多说，直接开始实验

这里我们采用的是单核CPU进行测试，型号为I7-9700K。需要用到的库只有numpy和scipy，当然主要的工具其实都在scipy的子库中。

接下来我们生成一个5000*5000的随机矩阵，并对它进行对称正定化，其结果记为矩阵大G。另外为了方便验证正确性，我们生成一个5000*1的向量。也就是说实际上我们是要解一个线性方程组。（这个做法是合理的，因为解线性方程组只比直接求逆多了一个矩阵乘向量的过程，其计算时间比求逆本身要小一个量级）。

另外可能有同学已经注意到上面的对称正定化操作有点眼熟悉。如果细心的同学应该会发现我们在讲牛顿迭代时其中一个变种的算法就用到了这种操作，这个算法叫做Levenberg-Marquart。

接下来开始实验。这里我们要用的算法分别是conjugate gradient, cholesky decomposition 以及LU分解。（scipy中的inv默认就是调用的LU分解）。

具体的计时方法非常简单，让程序运行100次，求出平均运行时间以及运行时间的STD。这样基本能保证运行结果相对稳定，对比也就更有力。图中我们只列出CG算法的代码，其它两种完全类似。

接下来看看结果：

首先是CG算法。平均运行时长为4000多毫秒，STD为32毫秒。总体来说十分稳定。

然后是cholesky。平均运行时长为555毫秒，STD不到5毫秒。非常优秀。

最后是lu分解。平均运行时长为1800多毫秒，STD80多秒。仍然十分稳定。（另外这里可能是因为我电脑其间断了一下网，导致它的STD出现了一点波动）。

这样一来就很容易看出优劣了。cholesky对于正定对称矩阵求逆的效率最高，LU次之，最慢是CG。

这还没完，再看看它们的精度如何：

精度对比非常简单，直接将它们解出的近似值带回线性方程组即可。这里可以看到cholesky和lu分解的误差都在10的负12次方量级上，非常精确。而CG则只小数点后一位，差距甚大。这里CG我们只用了它默认的精度tol=1e-6。如果再设小一些它的效率还会更低。那么综上来看，高下立判。

我们这里得出的结论和理论上的结论是完全一致的。LU分解是对于普通的矩阵通用的方法，因此其效率不是特别高。cholesky本身就是针对正定对称矩阵特有的方法，其效率最高。

可能许多同学会疑惑：CG通常在文献中被报导为解决大矩阵时最快的算法之一，为何在这里表现如此之差？这里我们不做过多的分析。简单提示一下：其实从算法的复杂度，再结合实际的情况，基本上可以猜到原因。

另外，还有一种情况我们也还没有测。那就是当矩阵特别大时，CG优势会不会明显一些？这个是有可能的。因为CG在整个运算过程中只存一个矩阵，而另外两个方法都需要存储额外的矩阵，即便对稀疏矩阵进行特殊处理，其复杂度也是难以显著降下来。不过这种测试太费时间和机器，有兴趣的同学可以自己去试一试。

本文源码都放在公从号 【mathit】上，关注后回复【正定对称矩阵】即可获取。