线性方程组的求解方法的选择和加速

Posted 陆嵩

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了线性方程组的求解方法的选择和加速相关的知识,希望对你有一定的参考价值。

线性方程组的求解方法的选择和加速

迭代法概述

求解线性方程组一般有直接解法(分解,如 LU,Cholesky 等)、迭代解法(Krylov、分块 Jacobi、多重网格等)和启发式方法。如果系数矩阵 A 很大并且是稀疏矩阵,分解方法一般情况下将不会有效。大部分情况下,小规模计算喜欢用直接解法,不存在时间问题了,就不需要加速。大规模问题的计算,我们青睐于用迭代解法,也有的求解器做出了适用于大规模问题的直接分解方法,比如 MUMPS, SuperLU 之类的,这是后话,我们主要谈迭代方法,直接解法和启发式解法就不说了。

迭代方法又包括矩阵分裂方法和 Krylov 子空间方法两大类。

矩阵分裂方法

矩阵分裂方法就是你们经常听到的 Jacobi 迭代、Gauss_Seidel 迭代和最佳因子 SOR 迭代等等。参见我的博文,

Jacobi迭代、Gauss_Seidel迭代和最佳因子SOR迭代的比较

SOR和SSOR迭代

共轭梯度方法(CG)MATLAB编程及其和Gauss_Seidel方法的一个比较

还有更多相关的文章,我就不列举了,你们可以在博客首页的搜索框里面搜索。对于一般的矩阵,这类算法不一定收敛,即使收敛,也有可能收敛得很慢。这个还是和问题有关。

Krylov 子空间方法

krylov 子空间方法就多了,简单列举几个,比如最速下降方法、共轭梯度方法、Lanczos方法、Arnoldi 方法、GMRES 方法(广义极小剩余法)、GCR 方法(广义共轭剩余法)、QR 方法(部分)、乘幂方法(子空间迭代)等等,不一而足。

具体的理论可以参考书籍:《矩阵计算的理论与方法》,徐树方,北京大学出版社。

GMRES中文叫广义极小剩余法,它是上世纪末到本世界初被人们认为是求解大型稀疏非对称线性方程组的最有效的方法之一。

对角占优的问题,我们可以生成质量更好的预条件子,适合预条件方法。

我这里也列出一小部分我的博文参考。

方程组求解的直接法与迭代法实现

共轭梯度(CG)方法

共轭梯度方法(CG)MATLAB编程及其和Gauss_Seidel方法的一个比较

最速下降法求特征值matlab程序

线性方程组迭代方法的选择

迭代方法可用于任何矩阵,但它们通常应用于直接求解速度慢的大型稀疏矩阵。鉴于很多人不想或者没能力进行编程实现迭代算法,我这里给出 MATLAB 中列出的一些可以调用的迭代求解函数,但是,MATLAB 的求解器非常具有局限性,仅供参考而已。

下面列出在 MATLAB 中最常见的一些求解方法,以及其适用的矩阵。

共轭梯度类

pcg(预条件共轭梯度法)

  • 系数矩阵必须是对称正定矩阵。
  • 由于只需要存储有限数量的向量,因此是最有效的对称正定方程组求解器。

bicg(双共轭梯度法)

  • 系数矩阵必须是方阵。
  • bicg 的计算成本低,但收敛不规则而且不可靠。
  • bicg 具有很重要的历史地位,因为很多迭代算法都是基于它改进而来。

bicgstab(双共轭梯度稳定法)

  • 系数矩阵必须是方阵。
  • 交替使用 BiCG 步骤与 GMRES(1) 步骤,以提高稳定性。

bicgstabl(双共轭梯度稳定法 (l))

  • 系数矩阵必须是方阵。
  • 交替使用 BiCG 步骤与 GMRES(2) 步骤,以提高稳定性。

cgs(共轭梯度二乘法)

  • 系数矩阵必须是方阵。
  • 需要与 bicg 相同的每次迭代操作数,但要避免通过操作平方残差来使用转置。

lsqr(最小二乘法)

  • 唯一适用于矩形方程组的求解器。
  • 在解析上等效于应用于标准方程 (A’*A)*x = A’*b 的共轭梯度法 (PCG)。

残差类

minres(最小残差法)

  • 系数矩阵必须是对称矩阵,但无需是正定矩阵。
  • 每次迭代都会最小化 2-范数中的残差,从而保证算法能够逐步取得进展。
  • 不会出现算法崩溃(算法无法逼近解而停止执行)。

symmlq(对称 LQ)

  • 系数矩阵必须是对称矩阵,但无需是正定矩阵。
  • 求解投影方程组,并使残差与所有先前的残差正交。
  • 不会出现算法崩溃(算法无法逼近解而停止执行)。

gmres(广义最小残差法)

  • 系数矩阵必须是方阵。
  • 最可靠的算法之一,因为每次迭代都能最小化残差范数。
  • 工作量和所需的存储量随迭代次数线性增加。
  • 选择适当的 restart 值至关重要,以避免不必要的工作和存储。

qmr(拟最小残差法)

  • 系数矩阵必须是方阵。
  • 每次迭代的开销略高于 bicg,但提供了更好的稳定性。

tfqmr(无转置拟最小残差法)

  • 系数矩阵必须是方阵。
  • 当内存有限时,是尝试用于求解对称不定方程组的最佳求解器。

更牛逼的加速手段:多重网格和并行(并发)计算(包含 MPI 加速,GPU 加速,多线程、多进程、共享内存譬如 OpenMP 等等)

多重网格方法的加速效果是非常客观的,迭代很少的步数,就可以快速地消除高频误差。可以参考我的博文,比如:

多重网格方法

并行和高并发的手段,可以采用多进程和多线程的方法,这可能需要对串行代码做一个并行化的设计,使得一些无关的操作可以在不同的节点或者线程下进行同时计算。这个就不多说了,暂略。

其他加速方法

其他的加速方法也有很多,就比较零散了,比如方程组求解的切比雪夫半迭代加速方法,可以参考我的博文。

方程组求解的切比雪夫半迭代加速方法

另外,非常值得强调的是,求解线性方程组,本质上求一个椭圆问题的极小化问题。那么,就可以用很多优化的手段来解决这个问题。比如采用拟牛顿方法、信赖域方法等等。我也给出几个示例的参考博文。

求解非约束优化问题的拟牛顿方法(BFGS、DFP)

求极小值的线搜索方法应用(SD,Newton)

MATLAB 求函数极值的内置函数一览表(实则优化算法函数汇总)

信赖域狗腿(dogleg)方法

强烈推荐大家仔细读一读我的两篇关于优化的类科普的基础知识文章。使用优化手段求解线性方程组,也是非常 nice 的。

数据科学和机器学习中的优化理论与算法(上)

数据科学和机器学习中的优化理论与算法(下)

以上是关于线性方程组的求解方法的选择和加速的主要内容,如果未能解决你的问题,请参考以下文章

线性方程理论说明和Eigen解线性方程求解方法汇总

求解线性代数矩阵方程

matlab中如何求解齐次线性方程组(代数矩阵)的非零解

SciPy 非线性方程求解 | Python技能树征题

Matlab求解线性方程组Ax=b的几种常见方法Matlab求解线性方程组Ax=b的几种常见方法

线性代数,求解矩阵方程