共轭梯度法（Conjugate Gradient）_共轭梯度法伪代码-CSDN博客

本文链接：https://blog.csdn.net/maum61/article/details/88088279

共轭梯度法（Conjugate Gradient）是介于最速下降法与牛顿法之间的一个方法，它仅需利用一阶导数信息，但克服了最速下降法收敛慢的缺点，又避免了牛顿法需要存储和计算Hesse矩阵并求逆的缺点，共轭梯度法不仅是解决大型线性方程组最有用的方法之一，也是解大型非线性最优化最有效的算法之一。

以上参考百度百科的共轭梯度法，但是一开始没看明白，结合《数值分析》中的一些解释，结合自己的理解，算是大概弄明白了。

wikipedia关于共轭梯度法描述的比较清晰。

一什么是A共轭？

令A是对称正定的n*n矩阵，对于n维的列向量v和w，定义A内积为：

$(v,w)_{A}=v^{T}Aw$

当 $(v,w)_{A}=0$ 时，称列向量v和w为A共轭。

特别的，当 $v=w\neq0$ 时， $(v,w)_{A}>0$

二共轭梯度法的伪代码分析

共轭梯度直接法主要作用是求解Ax=b的精确解，但是如果b的维数较高，直接法会消耗大量时间。而共轭梯度法作为迭代方法的话，对于求解高维问题的逼近解或稀疏解是很有用的，我们可以用部分的梯度方向去获得一个符合要求的逼近解，这样可以缓解直接共轭梯度方法的耗时问题。

我们构造 $\bigtriangledown f=Ax-b$ ，问题是求解Ax=b。也就是令 $\bigtriangledown f=Ax-b$ 为0，那么对 $\bigtriangledown f=Ax-b$ 积分，可得到二次泛函：

那问题也就转向了求f(x)的极小值点问题。我们知道，当多元函数的二次倒数也就是黑塞矩阵是正定的时候，那么驻点就是极小值。因此要求： $\bigtriangledown^{2} f=A$ 为正定的。

为了求f(x)的最小值，只要我们一步一步得使 $\bigtriangledown f=Ax-b$ 趋近于0即可。假设x的值是x0，那么梯度的误差为： $Ax_{0}-b$ ，因此我们要往反方向缩小误差。取初始的方向为 $\boldsymbol{p}_{0}=\boldsymbol{b}-A\mathbf{x}_{0}$ ，同时，p0也是残差。

令第k步的残差为 $\boldsymbol{r}_{k}=\boldsymbol{b}-A\mathbf{x}_{k}$

如何找函数f(x)下降最快的方向：

这里其实还可以这样理解，泛函在xk处的最优下降方向，不仅仅与负梯度有关，还应该和刚刚使用过的方向pk-1有关，那和二者的关系有多大呢？用 $\xi \eta$ 表示。所以对二者分别偏导，计算出来二者的值。

向量相减，得到的是与后者垂直的差向量。

Apk在计算步长的时候已经计算过了，因此可以将其保存，用于rk+1的计算。

伪代码如下：

$x_{0}=initial value$ %这是给临时解的一个初始值

$d_{0}=r_{0}=b-Ax_{0}$ %二次泛函的负梯度反方向

$x_{1}=x_{0}+\alpha _{0}d_{0}$ %计算沿负梯度方向的临时解

$r_{1}=b-Ax_{1}$ %这是利用初始方向计算的负梯度

for i=1:k %共轭梯度的新的搜索方向不再取负梯度。

$if r_{i}=0,stop,end$ %这里残差=0是有可能的，因为A的正定对称阵，因此是满秩矩阵，所以b必然在A的张成空间中

$\alpha_{i}=\frac{r_{i}^{T}r_{i}}{d_{i}^{T}Ad_{i}}$ %注意

$x_{i+1}=x_{i}+\alpha _{i}d_{i}$ %这里是更新临时解

$r_{i+1}=r_{i}-\alpha _{i}Ad_{i}$ %更新残差，这个是负梯度方向，

$\beta _{i}=\frac{r_{i+1}^{T}r_{i+1}}{r_{i}^{T}r_{i}}$ %选择系数 $\beta_{i}$

$d_{i+1}=r_{i+1}+\beta _{i}d_{i}$ %更新方向，共轭梯度方向使用当前负梯度方向和前一次方向的线性组合

end

该算法名叫“共轭梯度”，那就必然用到了共轭的概念。

首先，我们要明确，算法的目的，是要求解Ax=b，是一步一步迭代逼近精确解。

$\alpha_{i}=\frac{r_{i}^{T}r_{i}}{d_{i}^{T}Ad_{i}}$ 是如何得到的？首先，我们希望使用前一次的方向 $d_{i}$ 来更新本次的临时解： $x_{i+1}=x_{i}+\alpha _{i}d_{i}$ ，对其做变换：

$Ax_{i+1}=A(x_{i}+\alpha _{i}d_{i})\Rightarrow b-Ax_{i+1}=b-A(x_{i}+\alpha _{i}d_{i}) \Rightarrow r_{i+1}=r_{i}-\alpha _{i}Ad_{i}$

因为是循环迭代，并且我们规定，新的余项 $r_{i+1}$ 要和前面的所有余项都正交，这样能保证在n步以内彻底消除余项。那么我们要精心选择 $\alpha_{i}$ ，使得新的余项 $r_{i+1}$ 和方向 $d_{i}$ 正交。因此，前面得到的式子左乘 $d_{i}$ 得：

$d^{T}_{i}r_{i+1}=0=d^{T}_{i}r_{i}-\alpha_{i}Ad_{i}$

故 $\alpha_{i}=\frac{d_{i}^{T}r_{i}}{d_{i}^{T}Ad_{i}}$ ，这里得到的 $\alpha_{i}$ 与 $\alpha_{i}=\frac{r_{i}^{T}r_{i}}{d_{i}^{T}Ad_{i}}$ 不同，但是可以稍作转换得到。

新的余项 $r_{i+1}$ 要和前面的所有余项都正交，新的方向由该式得到：

$d_{i+1}=r_{i+1}+\beta _{i}d_{i}$

则有：

$d_{i}-r_{i}=\beta _{i-1}d_{i-1}$ 。因为新的余项 $r_{i+1}$ 要和前面的所有余项都正交，所以 $r_{i}\perp d_{i-1}\Rightarrow r_{i}\perp \beta_{i} d_{i-1}$

，此处 $\beta_{i}$ 是标量。所 $d_{i}-r_{i}=\beta _{i-1}d_{i-1}$ 左乘 $r^{T}_{i}$ 有：

$r^{T}_{i}\beta _{i-1}d_{i-1}=0=r^{T}_{i}d_{i}-r^{T}_{i}r_{i}\Rightarrow r^{T}_{i}d_{i}=r^{T}_{i}r_{i}\Rightarrow d^{T}_{i}r_{i}=r^{T}_{i}r_{i}$

所以得 $\alpha_{i}=\frac{d_{i}^{T}r_{i}}{d_{i}^{T}Ad_{i}}=\frac{r_{i}^{T}r_{i}}{d_{i}^{T}Ad_{i}}$

再说说如何更新 $\beta _{i}$

算法需要 $d_{i+1}=r_{i+1}+\beta _{i}d_{i}$ ，也就是新的方向由新的余项加上 $\beta _{i}$ 倍的原来的方向。

此处要保证所有的方向都是两两A共轭，也就是 $d^{T}_{i+1}Ad_{i}=0,d^{T}_{i}Ad_{i}>0$

所以对 $d_{i+1}=r_{i+1}+\beta _{i}d_{i}$ 构造共轭：

$d^{T}_{i}Ad_{i+1}=d^{T}_{i}A(r_{i+1}+\beta _{i}d_{i})=0\Rightarrow \beta _{i}=-\frac{d^{T}_{i}Ar_{i+1}}{d^{T}_{i}Ad_{i}}$

至于如何从 $\beta _{i}=-\frac{d^{T}_{i}Ar_{i+1}}{d^{T}_{i}Ad_{i}}$ 变到 $\beta _{i}=\frac{r_{i+1}^{T}r_{i+1}}{r_{i}^{T}r_{i}}$ ，书上也没提，我也实在参不透，希望知道的朋友留言相告，不吝赐教，谢谢！