Reference

Broyden–Fletcher–Goldfarb–Shanno algorithm

Sherman–Morrison formula

Woodbury matrix identitty

牛顿法与拟牛顿法学习笔记（四）BFGS 算法

Broyden类算法：BFGS算法的迭代公式推导（应用两次Sherman-Morrison公式）

Sherman-Morrison公式在BFGS算法的应用

BFGS算法中Sherman-Morrison-Woodbury公式的使用

【概述】

BFCS 算法是建立在阻尼牛顿法之上的，其以发明者 Broyden、Fletcher、Goldfarb、Shanno 四人姓名的首字母命名的，与 DFP 算法相比，其性能更佳，目前已成为求解无约束非线性优化问题的常用方法

对于阻尼牛顿法的搜索方向 $\mathbf{d_k}=-H_k^{-1}\cdot \mathbf{g_k}$，根据拟牛顿条件，BFCS 选用近似矩阵 $B_k$ 作为海森矩阵 $H_k$ 的近似，其迭代格式为：

$B_{k+1}=B_k+\triangle B_k,\quad k=0,1,2,...$

其中，$B_0$ 取单位矩阵 $I$，因此，算法的关键就是每一次迭代的校正矩阵 $\triangle B_k$ 如何构造

关于拟牛顿条件详见：拟牛顿迭代法

【BFCS 算法】

算法原理

假设校正矩阵 $\triangle B_k$ 与拟牛顿条件 $\mathbf{y_k} = B_{k+1}\boldsymbol{\delta_k}$、近似矩阵 $B_k$ 有关

采用待定系数法，设 $\alpha$、$\beta$ 为待定系数，向量 $\boldsymbol{\mu},\boldsymbol{\nu} \in \mathbb{R}^n$ 为待定列向量，则校正矩阵 $\triangle B_k$ 可待定为：

$\triangle B_k=\alpha\boldsymbol{\mu}\boldsymbol{\mu}^T+\beta\boldsymbol{\nu}\boldsymbol{\nu}^T$

从形式上看，这种待定公式保证了校正矩阵 $\triangle B_k$ 的对称性

结合拟牛顿条件 $\mathbf{y_k} = B_{k+1}\boldsymbol{\delta_k}$，将 $\triangle B_k$ 的待定式带入迭代式，有：

$\begin{align*} \mathbf{y_k} &= (B_k+\triangle B_k)\boldsymbol{\delta_k}\\ &= B_k\boldsymbol{\delta_k}+\alpha\boldsymbol{\mu}\boldsymbol{\mu}^T\boldsymbol{\delta_k}+\beta\boldsymbol{\nu}\boldsymbol{\nu}^T\boldsymbol{\delta_k} \notag \\ &= B_k\boldsymbol{\delta_k}+\boldsymbol{\mu}(\alpha\boldsymbol{\mu}^T\boldsymbol{\delta_k})+\boldsymbol{\nu}(\beta\boldsymbol{\nu}^T\boldsymbol{\delta_k}) \notag \\ &= B_k\boldsymbol{\delta_k}+(\alpha\boldsymbol{\mu}^T\boldsymbol{\delta_k})\boldsymbol{\mu}+(\beta\boldsymbol{\nu}^T\boldsymbol{\delta_k})\boldsymbol{\nu} \notag \end{align*}$

可以发现，括号中的 $\alpha\boldsymbol{\mu}^T\boldsymbol{\delta_k}$ 与 $\beta\boldsymbol{\nu}^T\boldsymbol{\delta_k}$ 是两个数，取最简单的情况，进行如下赋值：

$\left\{\begin{array}{rl} \alpha\boldsymbol{\mu}^T\boldsymbol{\delta_k} &=& 1 \\ \beta\boldsymbol{\nu}^T\boldsymbol{\delta_k} &=& -1 \end{array} \right.$

并将其带回 $\mathbf{y_k}$ 中，有：

$\mathbf{y_k}-B_k\boldsymbol{\delta_k}=\boldsymbol{\mu}-\boldsymbol{\nu}$

为使上式成立，直接取最简单情况 $\boldsymbol{\mu}=\boldsymbol{\delta_k}$，$\boldsymbol{\nu}=B_k\boldsymbol{\delta_k}$，此时，对于待定系数 $\alpha$、$\beta$ 有：

$\left\{\begin{array}{rl} \alpha &=& \frac{1}{\mathbf{y_k}^T\boldsymbol{\delta_k}} \\ \beta &=& -\frac{1}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}} \end{array} \right.$

这时，校正矩阵 $\triangle B_k$ 已构造出来，即：

$\triangle B_k=\frac{\mathbf{y_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\mathbf{\delta_k}}-\frac{B_k\boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB_k^T}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}}$

由于近似矩阵 $D_k$ 为近似成海森矩阵 $H_k$ 的正定矩阵，故有：

$B_k^T=B_k$

因此，BFCS 最终的迭代式为：

$B_{k+1}=B_k+\frac{\mathbf{y_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}}-\frac{B_k\boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB_k}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}},\quad k=0,1,2,...$

可以发现，对于 BFGS 与 DFP 相比，除了将近似矩阵 $D$ 换为 $B$ 外，其他地方只是将 $\boldsymbol{\delta_k}$ 替换为 $\mathbf{y_k}$，将 $\mathbf{y_k}$ 替换为 $\boldsymbol{\delta_k}$

算法流程

BFGS 算法的算法流程如下：

给定初值 $\mathbf{x_0}$ 和精度阈值 $\varepsilon$，并令 $B_0=I$，$k=0$
计算梯度向量 $\mathbf{g_k}$，与搜索方向 $\mathbf{d_k}=-B_k^{-1}\cdot \mathbf{g_k}$
利用 $\lambda_k=\arg \min\limits_{\lambda\in \mathbb{R}} f(\mathbf{x_k}+\lambda \mathbf{d_k})$ 得到步长 $\lambda_k$，并令 $\boldsymbol{\delta_k}=\lambda\mathbf{d_k}$，$\mathbf{x_{k+1}}=\mathbf{x_k} + \mathbf{\delta_k}$
若 $\mathbf{x_{k+1}}$ 的梯度向量 $\mathbf{g_{k+1}}$ 的 L1 范数 $||\mathbf{g_{k+1}}||$ 满足 $||\mathbf{g_{k+1}}||<\varepsilon$，则停止迭代，得近似解 $\mathbb{x}^*=\mathbf{x_{k+1}}$，否则转至步骤 5
计算两次迭代的梯度差 $\mathbf{y_k}=\mathbf{g_{k+1}}-\mathbf{g_k}$
计算第 $k+1$ 次的近似矩阵 $B_{k+1}=B_k+\frac{\mathbf{y_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}}-\frac{B_k\boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB_k}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}}$
令 $k=k+1$，转至步骤 2

【谢尔曼-莫里森公式下的 BFGS 算法】

算法原理

对于阻尼牛顿法来说，其搜索方向为：

$\mathbf{d_k}=-H_k^{-1}\cdot \mathbf{g_k}$

BFGS 算法选择使用近似矩阵 $B_k$ 来作为海森矩阵 $H_k$ 的近似，因此，这不可避免的需要在计算出 $B_k$ 后再求 $B_k^{-1}$，通常情况下，对于 $B_k^{-1}$ 的计算，都是通过解线性方程组 $B_k\mathbf{d_k}=-\mathbf{g_k}$ 进行的，这无疑影响了算法的运行速度

为此，在 BFGS 算法的基础上，使用谢尔曼-莫里森公式（Sherman-Morrison Formula），将 BFGS 的迭代式进行改进，直接给出 $B_{k+1}^{-1}$ 与 $B_k^{-1}$ 的关系式，即：

$B_{k+1}^{-1} = (I-\frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}})B_k^{-1}(I-\frac{\mathbf{y_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}})+\frac{\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}}$

从而对 BFGS 算法进一步优化，使得整个算法不再需要解线性方程组，仅使用矩阵-向量运算即可完成

为算法流程的清晰性，不再出现矩阵求逆符号，统一将 $B_k^{-1}$ 换为 $G_k$，并称 $G_k$ 为特殊 BFGS 矩阵，即：

$G_{k+1} =(I-\frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}})G_k(I-\frac{\mathbf{y_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}})+\frac{\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}}$

算法流程

此时 BFGS 算法流程如下：

给定初值 $\mathbf{x_0}$ 和精度阈值 $\varepsilon$，并令 $G_0=I$，$k=0$
计算梯度向量 $\mathbf{g_k}$，与搜索方向 $\mathbf{d_k}=-G_k\cdot \mathbf{g_k}$
利用 $\lambda_k=\arg \min\limits_{\lambda\in \mathbb{R}} f(\mathbf{x_k}+\lambda \mathbf{d_k})$ 得到步长 $\lambda_k$，并令 $\boldsymbol{\delta_k}=\lambda\mathbf{d_k}$，$\mathbf{x_{k+1}}:=\mathbf{x_k} + \boldsymbol{\delta_k}$
若 $\mathbf{x_{k+1}}$ 的梯度向量 $\mathbf{g_{k+1}}$ 的 L1 范数 $||\mathbf{g_{k+1}}||$ 满足 $||\mathbf{g_{k+1}}||<\varepsilon$，则停止迭代，否则转至步骤 5
计算两次迭代的梯度差 $\mathbf{y_k}=\mathbf{g_{k+1}}-\mathbf{g_k}$
计算第 $k+1$ 次的近似矩阵 $G_{k+1} = (I-\frac{\mathbf{\delta_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\mathbf{\delta_k}})G_k(I-\frac{\mathbf{y_k}\mathbf{\delta_k}^T}{\mathbf{y_k}^T\mathbf{\delta_k}})+\frac{\mathbf{\delta_k}\mathbf{\delta_k}^T}{\mathbf{y_k}^T\mathbf{\delta_k}}$
令 $k=k+1$，转至步骤 2

算法推导

谢尔曼-莫里森公式

谢尔曼-莫里森公式（Sherman-Morrison Formula）是以 Sherman 和 Morrison 命名，在线性代数中，是求解逆矩阵的一种方法，其内容如下：

设 $A\in \mathbb{R}^{n\times n}$ 是非奇异方阵，$\boldsymbol{\mu},\boldsymbol{\nu} \in \mathbb{R}^{n}$，则 $(A+\boldsymbol{\mu}\boldsymbol{\nu}^T)$ 可逆，当且仅当 $1+\boldsymbol{\nu}^TA\boldsymbol{\mu}\neq0$，其逆矩阵为：
$(A+\boldsymbol{\mu}\boldsymbol{\nu}^T)^{-1}=A^{-1}-\frac{A^{-1}\boldsymbol{\mu}\boldsymbol{\nu}^TA^{-1}}{1+\boldsymbol{\nu}^TA^{-1}\boldsymbol{\mu}}$

通过使用两次 Sherman-Morrison 公式，可以在 BFGS 算法迭代式的基础上得到蕴含 $B_{k+1}^{-1}$ 与 $B_k^{-1}$ 的关系式的 BFCS 算法的迭代式，证明如下：

已知 BFGS 算法迭代式为：

其中，近似矩阵 $B_k\in \mathbb{R}^{n\times n}$ 是可逆的实对称矩阵，梯度差 $\mathbf{y_k}=\mathbf{g_{k+1}}-\mathbf{g_k}\in \mathbb{R}^n$，$\mathbf{x}$ 的增量 $\boldsymbol{\delta_k}=\mathbf{x_{k+1}}-\mathbf{x_k} \in \mathbb{R}^n$

对于 $-\frac{B_k\boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB_k}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}}$，有：

$-\frac{B_k\boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB_k}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}} =\frac{1}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}}(B_k\boldsymbol{\delta_k})(B_k\boldsymbol{\delta_k})^T$

其中，$\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}$ 为二次型形式的标量，$B_k\boldsymbol{\delta_k}\in \mathbb{R}^n$

记 $A=B_{k}+\frac{\mathbf{y_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}}$，那么有：

$B_{k+1}^{-1}=A+\frac{1}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}}(B_k\boldsymbol{\delta_k})(B_k\boldsymbol{\delta_k})^T$

根据 Sherman-Morrison 公式，有：

$\begin{align} B_{k+1}^{-1} &= \Big(A+\frac{1}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}}(B_k\boldsymbol{\delta_k})(B_k\boldsymbol{\delta_k})^T \Big)^{-1} \notag \\ &= A^{-1} - \frac{A^{-1} \Big( -\frac{1}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}} (B_k\boldsymbol{\delta_k}) (B_k\boldsymbol{\delta_k})^T \Big) A^{-1}}{1+ \Big(-\frac{1}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}} (B_k\boldsymbol{\delta_k})^T A^{-1} (B_k\boldsymbol{\delta_k})\Big)} \notag \\ &= A^{-1}+A^{-1}\frac{B_k\boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB_k}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k} - \boldsymbol{\delta_k}^TB_kA^{-1}B_k\boldsymbol{\delta_k}}A^{-1} \notag \end{align}$

对于 $A=B_{k}+\frac{\mathbf{y_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}}$，分母 $\mathbf{y_k}^T\boldsymbol{\delta_k}$ 为标量，再次利用 Sherman-Morrison 公式，有：

$\begin{align} A^{-1} &= \Big( B_k+\frac{\mathbf{y_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}} \Big)^{-1} \notag \\ &= B_k^{-1}-\frac{B_k^{-1} \frac{\mathbf{y_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}} B_k^{-1}}{1+ \frac{1}{\mathbf{y_k}^T\boldsymbol{\delta_k}}\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}} \notag \\ &= B_k^{-1}-\frac{B_k^{-1}\mathbf{y_k}\mathbf{y_k}^TB_k^{-1}}{\mathbf{y_k}^T\boldsymbol{\delta_k}+\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}} \notag \end{align}$

其中，分母 $\mathbf{y_k}^T\boldsymbol{\delta_k}$ 是内积形式的标量，$\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}$ 是二次型形式的标量，记：

$Q=\mathbf{y_k}^T\boldsymbol{\delta_k}+\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}$

则有：

$A^{-1}=B_k^{-1}-\frac{B_k^{-1}\mathbf{y_k}\mathbf{y_k}^TB_k^{-1}}{Q}$

对于 $B_{k+1}^{-1}= A^{-1}+A^{-1}\frac{B_k\boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB_k}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k} - \boldsymbol{\delta_k}^TB_kA^{-1}B_k\boldsymbol{\delta_k}}A^{-1} $，先对后半部分进行化简，有：

$\begin{align} 原式 &= A^{-1}\frac{B_k\boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB_k}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k} - \boldsymbol{\delta_k}^TB_k \Big(B_k^{-1} -\frac{B_k^{-1}\mathbf{y_k}\mathbf{y_k}^TB_k^{-1}}{Q}\Big) B_k\boldsymbol{\delta_k}}A^{-1} \notag \\ &= A^{-1} \frac{B_k\boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB_k}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k} - \boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}+ \frac{\boldsymbol{\delta_k}^T\mathbf{y_k}\mathbf{y_k}^T\boldsymbol{\delta_k}}{Q}} A^{-1} \notag \\ &= \frac{Q}{\boldsymbol{\delta_k}^T\mathbf{y_k}\mathbf{y_k}^T\boldsymbol{\delta_k}} \Big( B_k^{-1}-\frac{B_k^{-1}\mathbf{y_k}\mathbf{y_k}^TB_k^{-1}}{Q} \Big) B_k \boldsymbol{\delta_k}\boldsymbol{\delta_k}^T B_k A^{-1} \notag \\ &= \frac{Q}{\boldsymbol{\delta_k}^T\mathbf{y_k}\mathbf{y_k}^T\boldsymbol{\delta_k}} \Big( \boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB -\frac{B_k^{-1}\mathbf{y_k}\mathbf{y_k}^T \boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB_k}{Q} \Big) A^{-1} \notag \\ &\xlongequal[m=\mathbf{y_k}^T\boldsymbol{\delta_k}]{m=\boldsymbol{\delta_k}^T\mathbf{y_k}} \frac{Q}{m^2} \Big( \boldsymbol{\delta_k}\boldsymbol{\delta_k}^T - \frac{m\boldsymbol{\delta_k}\mathbf{y_k}^TB_k^{-1}}{Q} - \frac{B_k^{-1}m\mathbf{y_k}\boldsymbol{\delta_k}^T}{Q} + \frac{B_k^{-1}m^2\mathbf{y_k}\mathbf{y_k}^TB_k^{-1}}{Q^2} \Big) \notag \\ &= \frac{Q\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{m^2} - \frac{B_k^{-1}\boldsymbol{\delta_k}\mathbf{y_k}^T}{m} - \frac{B_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}^T}{m} + \frac{B_k^{-1}\mathbf{y_k}\mathbf{y_k}^TB_k^{-1}}{Q} \notag \end{align}$

此时，有：

$B_{k+1}^{-1} = A^{-1} + \frac{Q\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{m^2} - \frac{B_k^{-1}\boldsymbol{\delta_k}\mathbf{y_k}^T}{m} - \frac{B_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}^T}{m} + \frac{B_k^{-1}\mathbf{y_k}\mathbf{y_k}^TB_k^{-1}}{Q}$

再次将 $A^{-1}=B_k^{-1}-\frac{B_k^{-1}\mathbf{y_k}\mathbf{y_k}^TB_k^{-1}}{Q}$ 代入，化简可得：

$\begin{align} B_{k+1}^{-1} &= B_{k}^{-1} + \frac{Q\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{m^2} - \frac{\boldsymbol{\delta_k}\mathbf{y_k}^TB_k^{-1}}{m} - \frac{B_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}^T}{m} \notag \\ &= \Big( I-\frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{m}\Big)B_k^{-1} - \frac{B_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}^T}{m} + \frac{(\mathbf{y_k}^T\boldsymbol{\delta_k}+\mathbf{y_k}^TB_k^{-1}\mathbf{y_k})\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{m^2} \notag \\ &= \Big( I-\frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{m}\Big)B_k^{-1} - \frac{B_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}^T}{m} + \frac{\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{m} + \frac{\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{m^2} \notag \end{align}$

对于最后一项的分子 $\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T$，由于 $\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}$ 是二次型形式的标量，因此有：

$\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T = \boldsymbol{\delta_k}\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}^T$

故有：

$\begin{align} B_{k+1}^{-1} &= \Big( I-\frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{m}\Big)B_k^{-1} - \frac{B_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}^T}{m} + \frac{\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{m} + \frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{m} \frac{B_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}^T}{m} \notag \\ &= \Big( I-\frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{m}\Big)B_k^{-1} - \Big( I- \frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{m} \Big) B_k^{-1}\frac{\mathbf{y_k}\boldsymbol{\delta_k}^T}{m} + \frac{\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{m} \notag \\ &= \Big( I-\frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{m}\Big)B_k^{-1}\Big( I- \frac{\mathbf{y_k}\boldsymbol{\delta_k}^T}{m} \Big) +\frac{\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{m} \notag \end{align}$

将 $m=\mathbf{y_k}^T\boldsymbol{\delta_k}$ 带回，有：

$B_{k+1}^{-1}=\Big( I-\frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}}\Big) B_k^{-1} \Big( I- \frac{\mathbf{y_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}} \Big) +\frac{\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}}$

再记 $G_k=B_k^{-1}$，则有：

$G_{k+1}=\Big( I-\frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}}\Big) B_k^{-1} \Big( I- \frac{\mathbf{y_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}} \Big) +\frac{\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}}$

证明完毕

伍德伯里矩阵恒等式

伍德伯里矩阵恒等式（Woodbury Matrix Identity）是 Woodbury 在 Sherman-Morrison 公式的基础上进行的一般化推广，推广后的 Sherman-Morrison 公式一般称为 Sherman-Morrison-Woodbury 公式，其内容如下：

设 $A\in \mathbb{R}^{n\times n}$ 是非奇异方阵，$U,V \in \mathbb{R}^{n\times p}$，其中 $1\leq p \leq n$，则矩阵 $A+UV^T$ 可逆，当且仅当矩阵 $I+V^TA^{-1}U$ 可逆，且有：
$(A+UV^T)^{-1}=A^{-1}-A^{-1}U(I+V^TA^{-1}U)^{-1}V^TA^{-1}$

对于 BFCS 算法，通过使用两次 Sherman-Morrison 公式，可以在 BFGS 算法迭代式的基础上得到蕴含 $B_{k+1}^{-1}$ 与 $B_k^{-1}$ 的关系式的 BFCS 算法的迭代式

当使用 Sherman-Morrison-Woodbury 公式时，则仅需使用一次即可得到蕴含 $B_{k+1}^{-1}$ 与 $B_k^{-1}$ 的关系式的 BFCS 算法的迭代式，证明如下：

已知 BFGS 算法迭代式为：

$B_{k+1} = B_k + \frac{\mathbf{y_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}}-\frac{B_k\boldsymbol{\delta_k}\boldsymbol{\delta_k}^TB_k}{\boldsymbol{\delta_k}^TB_k\boldsymbol{\delta_k}},\quad k=0,1,2,...$

将后两项写成矩阵相乘的形式，有：

$B_{k+1} = B_k + \begin{bmatrix} -\frac{B_k\boldsymbol{\delta_k}}{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}} & \frac{\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}} \end{bmatrix} \begin{bmatrix} \boldsymbol{\delta_k}^T B_k \\ \mathbf{y_k}^T \end{bmatrix}$

记 $U=\begin{bmatrix}-\frac{B_k\boldsymbol{\delta_k}}{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}} & \frac{\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}}\end{bmatrix}$，$V=\begin{bmatrix}B_k\boldsymbol{\delta_k} &\mathbf{y_k}^T\end{bmatrix}$，则有：

$B_{k+1}=B_k+UV^T$

根据 Sherman-Morrison-Woodbury 公式，可得：

$B_{k+1}^{-1} = B_k^{-1} - B_k^{-1}U (I+V^TB_k^{-1}U)^{-1} V^TB_k^{-1}$

对于 $I+V^TB_k^{-1}U$，有：

$\begin{align} I+V^TB_k^{-1}U &= I + \begin{bmatrix} \boldsymbol{\delta_k}^T B_k \\ \mathbf{y_k}^T \end{bmatrix} B_k^{-1} \begin{bmatrix} -\frac{B_k\boldsymbol{\delta_k}}{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}} & \frac{\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}} \end{bmatrix} \notag \\ &= I + \begin{bmatrix} \boldsymbol{\delta_k}^T \\ \mathbf{y_k}^T B_k^{-1} \end{bmatrix} \begin{bmatrix} -\frac{B_k\boldsymbol{\delta_k}}{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}} & \frac{\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}} \end{bmatrix} \notag \\ &= I + \begin{bmatrix} -1 & 1 \\ -\frac{\mathbf{y_k}^T\boldsymbol{\delta_k}}{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}} & \frac{\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}} \end{bmatrix}\notag \\ &= \begin{bmatrix} 0 & 1 \\ -\frac{\mathbf{y_k}^T\boldsymbol{\delta_k}}{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}} & \frac{\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}}+1 \end{bmatrix} \notag \end{align}$

对于二阶矩阵来说，若 $\begin{bmatrix}a&b\\c&d\end{bmatrix}$ 可逆，那么有：

$\begin{bmatrix}a&b\\c&d\end{bmatrix} = \frac{1}{ad-bc} \begin{bmatrix}d&-b\\-c&a\end{bmatrix}$

于是，对于 $I+V^TB_k^{-1}U$ 来说，就有：

$\begin{align} (I+V^TB_k^{-1}U)^{-1} &=\frac{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}}{\mathbf{y_k}^T\boldsymbol{\delta_k}} \begin{bmatrix} \Big(\frac{\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}}+1\Big) & -1\\ \frac{\mathbf{y_k}^T\boldsymbol{\delta_k}}{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}} & 0 \end{bmatrix} \notag \\ &= \begin{bmatrix} \frac{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}}{\mathbf{y_k}^T\boldsymbol{\delta_k}}\Big(\frac{\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}}+1\Big) & -\frac{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}}{\mathbf{y_k}^T\boldsymbol{\delta_k}}\\ 1 & 0 \end{bmatrix} \notag \end{align}$

那么，对于 $B_k^{-1}U (I+V^TB_k^{-1}U)^{-1} V^TB_k^{-1}$，有：

$\begin{align} 原式 &= B_k^{-1} \begin{bmatrix}-\frac{B_k\boldsymbol{\delta_k}}{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}} & \frac{\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}}\end{bmatrix} \begin{bmatrix} \frac{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}}{\mathbf{y_k}^T\boldsymbol{\delta_k}}\Big(\frac{\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}}+1\Big) & -\frac{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}}{\mathbf{y_k}^T\boldsymbol{\delta_k}}\\ 1 & 0 \end{bmatrix} \begin{bmatrix} \boldsymbol{\delta_k}^T B_k \\ \mathbf{y_k}^T \end{bmatrix} B_k^{-1} \notag \\ &= \begin{bmatrix}-\frac{\boldsymbol{\delta_k}}{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}} & \frac{B_k^{-1}\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}} \end{bmatrix} \begin{bmatrix} \frac{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}}{\mathbf{y_k}^T\boldsymbol{\delta_k}}\Big(\frac{\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}}+1\Big) & -\frac{\boldsymbol{\delta_k}^T B_k \boldsymbol{\delta_k}}{\mathbf{y_k}^T\boldsymbol{\delta_k}}\\ 1 & 0 \end{bmatrix} \begin{bmatrix} \boldsymbol{\delta_k}^T \\ \mathbf{y_k}^TB_k^{-1} \end{bmatrix} \notag \\ &= \begin{bmatrix} -\frac{\boldsymbol{\delta_k}\Big(\frac{\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}}+1\Big)}{\mathbf{y_k}^T\boldsymbol{\delta_k}}+\frac{B_k^{-1}\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}} & \frac{\boldsymbol{\delta_k}}{\mathbf{y_k}^T\boldsymbol{\delta_k}} \end{bmatrix} \begin{bmatrix} \boldsymbol{\delta_k}^T \\ \mathbf{y_k}^TB_k^{-1} \end{bmatrix} \notag \\ &= \frac{\Big(\frac{\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}}+1\Big)\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T} {\mathbf{y_k}^T\boldsymbol{\delta_k}} + \frac{B_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}} + \frac{\boldsymbol{\delta_k}\mathbf{y_k}^TB_k^{-1}}{\mathbf{y_k}^T\boldsymbol{\delta_k}} \notag \end{align}$

综上，对于原式 $B_{k+1}^{-1} = B_k^{-1} - B_k^{-1}U (I+V^TB_k^{-1}U)^{-1} V^TB_k^{-1}$，有：

$\begin{align} B_{k+1}^{-1} &= B_k^{-1}-\Big(\frac{\mathbf{y_k}^TB_k^{-1}\mathbf{y_k}}{\mathbf{y_k}^T \boldsymbol{\delta_k}}+1\Big)\frac{\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T} {\mathbf{y_k}^T\boldsymbol{\delta_k}} - \frac{B_k^{-1}\mathbf{y_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}} - \frac{\boldsymbol{\delta_k}\mathbf{y_k}^TB_k^{-1}}{\mathbf{y_k}^T\boldsymbol{\delta_k}} \notag \\ &= \Big( I - \frac{\boldsymbol{\delta_k}\mathbf{y_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}} \Big) B_k^{-1} \Big( I - \frac{\mathbf{y_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}} \Big) + \frac{\boldsymbol{\delta_k}\boldsymbol{\delta_k}^T}{\mathbf{y_k}^T\boldsymbol{\delta_k}} \notag \end{align}$

再记 $G_k=B_k^{-1}$，则有：

证明完毕