预处理共轭梯度算法(Preconditioned Conjugate Gradients Method)

预处理共轭梯度算法(Preconditioned Conjugate Gradients Method)

给出百度百科上的解释：

预处理共轭梯度法

预处理共轭梯度法是。不必预先估计参数等特点。
共轭梯度法近年来在求解大型稀疏方程组中取得了较好的成效。理论上普通的共扼梯度法对于对称超正定方程，只要迭代步数达到方程的阶数就可以得到精确解，但实际上当系数矩阵的条件数(最大最小特征值之比)很大时，普通的共轭梯度法收敛速度很慢。预处理共轭梯度法对系数矩阵作预处理，以加速迭代收敛速度。

这个预处理共轭梯度算法，适用的环境是“正定的大型稀疏矩阵”，并且系数矩阵的条件数(最大最小特征值之比)很大的情况。

一般在机器学习中，我目前接触的问题中其实并不太会用到这个预处理共轭梯度算法，标准的共轭梯度算法就足以处理大多数的问题了。

共轭梯度法，在之前的博客中已经多次介绍并给出了对应的计算代码，这里就不具体介绍了，这里只讲一下这个预处理。

共轭梯度法，就是求解方程：Ax=b
其中，A矩阵为正定矩阵。

而预处理共轭梯度法，则是对矩阵A进行一个预处理，因为如果A是一个比较大的稀疏矩阵，并且A的系数矩阵的条件数(最大最小特征值之比)很大，那么即使使用共轭梯度法也需要较长的运算时间，因此可以在这种情况下，可以通过对A矩阵进行一个预处理得到等价的B矩阵，即Bx=b，这里的x和Ax=b中的x相同。

给出百度文库上的资料：

预处理共轭梯度法：

PS. 预处理方法有对角线预处理，不完全Cholesky分解预处理等。其目标就是把矩阵转换为矩阵B，并保证Ax=b，Bx=b，并且两者的x相同。这个预处理共轭梯度法更多的是用在物理学领域，在信息学中的应用还是比较有限的，因此只需要做到了解即可，不用强求掌握。

本文并没有找到具体的“预处理共轭梯度”的代码，但是给出了下面相关的资料，所需要者可以根据下面的资料执行实现这个预处理的步骤：

• 预处理共轭梯度法（1）
• 预处理共轭梯度法（2）
• 预处理（Preconditioning）
• 共轭梯度法（二）-预处理共轭梯度法

---

附录：

共轭梯度法，代码：

def cg(f_Ax, b, cg_iters=10, callback=None, verbose=False, residual_tol=1e-10):
    """
    Demmel p 312
    """
    p = b.copy()
    r = b.copy()
    x = np.zeros_like(b)
    rdotr = r.dot(r)
 
    fmtstr = "%10i %10.3g %10.3g"
    titlestr = "%10s %10s %10s"
    if verbose: print(titlestr % ("iter", "residual norm", "soln norm"))
 
    for i in range(cg_iters):
        if callback is not None:
            callback(x)
        if verbose: print(fmtstr % (i, rdotr, np.linalg.norm(x)))
        z = f_Ax(p)
        v = rdotr / p.dot(z)
        x += v * p
        r -= v * z
        newrdotr = r.dot(r)
        mu = newrdotr / rdotr
        p = r + mu * p
 
        rdotr = newrdotr
        if rdotr < residual_tol:
            break
 
    if callback is not None:
        callback(x)
    if verbose: print(fmtstr % (i + 1, rdotr, np.linalg.norm(x)))  # pylint: disable=W0631
    return x