本文只考虑最简单的一维问题，包括求根和求最小值点两个子问题：

  • 非线性方程求根 \(f(x)=0\)
  • 非线性方程求最小值点 \(x_0 = \text{argmin} f(x)\)
    • 精确一维搜索
    • 非精确一维搜索

由于线性问题是平凡的，因此默认\(f(x)\)是非线性的。

非线性方程求根

对分法

假设\(f(x)\)在有限区间\([a,b]\)满足：\(f(a)f(b)<0\) 异号，那么\(f(x)\)在\((a,b)\)至少有一个零点 \(x^*\)。

算法：

  1. 起始的计算区间 \([x_l,x_r]=[a,b]\)，此时满足两端异号 \(f(x_l)f(x_r)<0\) 。
  2. 选取中点 \(x_m = \frac12(x_l+x_r)\) ，计算\(f(x_m)\)。
  3. 如果中点 \(|f(x_m)| \le \varepsilon\)，或者区间长度\(|x_r-x_l|\le \varepsilon\) 已经足够小，则认为中点已经足够精确，已经找到零点，求解完成。
  4. 否则左右两个区间 \([x_l,x_m]\) ，\([x_m,x_r]\) 恰有一个满足两端异号，选取它作为新的计算区间\([x_l',x_r']\)，回到第二步。

收敛性：进行 \(k\) 次对分之后，计算区间长度为 \(2^{-k}(b-a)\)，因此至多进行 \(N\) 次对分。 \[ \frac{b-a}{2^N} \le \varepsilon \Rightarrow N = \left[\log_2(\frac{b-a}{\varepsilon})\right] \]

误差：第 k 次对分之后中点 \(x_m^{(k)}\) 满足 \[ \Vert x_m^{(k)} - x^*\Vert \le \frac{b-a}{2^{k+1}} \]

Remark：算法的前提条件保证了至少有一个零点，但如果在区间中有多个零点，则可能收敛到其中任一零点。

牛顿法

假设 \(f(x)\) 在 \(x\) 附近有一个零点 \(x^*\) ，牛顿法的流程是：从 \(x^{(0)}=x\) 出发，构造迭代点列 \(\{x^{(k)}\}\)，希望极限趋于零点 \[ \lim_{k \to \infty} x^{(k)} = x^* \] 具体的构造方法为：在\((x^{(k)},f(x^{(k)}))\)处作 \(f\) 的切线，切线与 x轴交点记作 \(x^{(k+1)}\)。 \[ \begin{aligned} 0 - f(x^{(k)}) &= f'(x^{(k)})(x^{(k+1)}-x^{(k)})\\ x^{(k+1)} &= x^{(k)} - \frac{f(x^{(k)})}{f'(x^{(k)})} \end{aligned} \] 迭代的终止条件有两个： \[ \begin{aligned} |f(x^{(k)})| \le \varepsilon_1 \\ \Vert x^{(k)} - x^{(k-1)})\Vert \le \varepsilon_2 \end{aligned} \]

Remark:

  • 牛顿法对迭代的初值很敏感，要求初值与真正的零点足够靠近（靠近程度也与 \(f\) 自身性质有关）；否则可能跳转并收敛到另外的零点，或者直接不收敛。
  • 在初值合适选取时，并且 \(f \in C^2\) 在 \(x^*\) 有单根时，收敛速度可以达到二阶。（如果在 \(x^*\) 是重根，可能需要对格式进行修改）

\[ \Vert x^{(k+1)}-x^*\Vert \le C \Vert x^{(k)}-x^*\Vert^2 \]

  • 对于单调且凸的函数，才足够保证牛顿法的全局收敛性，例如 \(f(x)=x^2\) 对应的是经典的平方根算法。

收敛速度的证明： 记 \(e^{(k)} = x^{(k)} - x^*\)，则有 \[ \begin{aligned} e^{(k+1)} &= x^{(k+1)} - x^* = x^{(k)} - \frac{f(x^{(k)})}{f'(x^{(k)})} - x^* \\ &= e^{(k)} - \frac{f(x^{(k)})}{f'(x^{(k)})} = \frac{e^{(k)} f'(x^{(k)}) - f(x^{(k)})}{f'(x^{(k)})} \end{aligned} \] 对\(f(x^{*}) = f(x^{(k)}-e^{(k)})\)在 \(x^{(k)}\) 处进行泰勒展开 \[ 0 = f(x^*) = f(x^{(k)}) - f'(x^{(k)}) e^{(k)} + \frac12 f''(\xi) (e^{(k)})^2 \] 因此有 \[ e^{(k+1)} = \frac12 \frac{f''(\xi)}{f'(x^{(k)})} (e^{(k)})^2 \approx C (e^{(k)})^2 \]

例：经典的平方根算法：计算 \(f(x)=x^2-R\) 的零点 \[ x^{(k+1)} = x^{(k)} - \frac{(x^{(k)})^2-R}{2x^{(k)}} = \frac12 \left( x^{(k)} + \frac{R}{x^{(k)}}\right) \]

Remark： 关于非线性方程求根的牛顿法，可以直接推广到非线性方程组上，求 \(X\in R^n\)，满足 \(F(X) = 0 \in R^n\)，这里为了简化问题，取约束个数和自变量维数相同，并且假设问题的解存在。最简单是线性情形 \(F(X)=AX\)，\(A \in R^{n\times n}\) 是常系数方阵。

迭代公式为 \[ X^{(k+1)} = X^{(k)} - J(X^{(k)})^{-1}F(X^{(k)}) \] 其中 \(J(X)\) 为 Jacobi 矩阵。 \[ J(X)_{ij} = \frac{\partial F_i}{\partial x_j}(X) \] 这里要求 \(F(X)\) 足够光滑，并且假设 Jacobi 矩阵始终可逆。

割线法

牛顿法需要使用导函数信息，有时对于非线性函数 \(f(x)\) 直接求导很困难，可以使用最近两个点的割线斜率进行近似。

\[ \begin{aligned} x^{(k+1)} &= x^{(k)} - \frac{f(x^{(k)})(x^{(k)}-x^{(k-1)})}{f(x^{(k)})-f(x^{(k-1)})} \end{aligned} \] 收敛的速度会略有降低（阶数\(\frac{1+\sqrt{5}}2\approx 1.62\)），此时为了启动格式，我们需要提供两个初始值\(x_0,x_1\)。

不动点法

对于方程 \(f(x)=0\) 转化为如下形式 \[ x = F(x) \] 即 \(f\) 的零点 \(x\) 是\(F\) 的不动点。然后可以据此设计迭代格式 \[ x^{(k+1)} = F(x^{(k)}) \]

基于压缩映射原理，我们希望在区间\([a,b]\)上的 \(F\) 是压缩映射，要求\(F:[a,b]\to [a,b]\) ，并且存在 \(0<L<1\) 使得 \[ \Vert F(x) - F(y) \Vert \le L \Vert x-y \Vert \]

牛顿法实质是一种不动点法： \[ F_N(x) = x - \frac{f(x)}{f'(x)} \] 割线法实质也是一种不动点法： \[ F_G(x) = x - \frac{f(x)^2}{f(x+f(x))-f(x)} \]

不动点法提供了非常自由的构造方式：针对具体问题，设计不同的 \(F\) 并检验压缩映射条件，就可以给出很多不动点迭代格式。

精确一维搜索

精确一维搜索通常要求 \(f\) 在区间 \([a,b]\) 上为单峰的函数，已知存在极小值，否则会计算失败。

黄金分割法

这个方法只需要计算 \(f(x)\) 的值，不需要导数信息。 主要的步骤是不断压缩计算区间，假设初始计算区间 \([x_l,x_r]\)，区间长度 \(x_r-x_l\)，比例系数 \(\rho \in (0,\frac12)\) 待定。

取内部左右两侧的两个点 \(a_l,a_r\) 为 \[ \begin{aligned} a_l = x_l + \rho (x_r-x_l)\\ a_r = x_r - \rho (x_r-x_l)\\ \end{aligned} \]

  • 如果 \(f(a_l)<f(a_r)\)，那么极小点位于 \(a_r\) 左侧，选取新的计算区间为 \([x_l',x_r']=[x_l,a_r]\)。
  • 如果 \(f(a_l)\ge f(a_r)\)，那么极小点位于 \(a_l\) 右侧，选取新的计算区间为 \([x_l',x_r']=[a_l,x_r]\)。

新的计算区间长度为 \((1-\rho)(x_r-x_l)\)。

我们希望可以充分利用计算的点值，以第一种情况为例，希望 \(a_l\) 可以在新计算区间 \([x_l',x_r']\) 也被利用，作为右侧点，需要满足 \[ \frac{1-2\rho}{1-\rho} = \rho \] 解得 \(\rho = \frac{3-\sqrt5}2 \approx 0.382\)。

（黄金分割法的思路是固定 \(\rho\)， 另一种思路是动态调整 \(\rho\)，取一个与斐波那契数列相关的系数 \(\rho_k\)，得到斐波那契数列法）

每次计算区间的压缩比例为 \(1-\rho \approx 0.618\) 即黄金分割比。 除了第一次压缩区间之外，每一次压缩区间只需要选择一个内部点进行求值。

算法：

  1. 起始的计算区间 \([x_l,x_r]=[a,b]\)。
  2. 选取左侧点 \(a_l = x_l + \rho(x_r-x_l)\) ，右侧点 \(a_r = x_r - \rho(x_r-x_l)\)，计算 \(f(a_l)\) 和 \(f(a_r)\)。
  3. 压缩区间：
    1. 如果 \(f(a_l)<f(a_r)\)，那么极小点位于 \(a_r\) 左侧，选取新的计算区间为 \([x_l',x_r']=[x_l,a_r]\)。\(a_l\) 可以作为下一步的 \(a_r'\)。
    2. 如果 \(f(a_l)\ge f(a_r)\)，那么极小点位于 \(a_l\) 右侧，选取新的计算区间为 \([x_l',x_r']=[a_l,x_r]\)。\(a_r\) 可以作为下一步的 \(a_l'\)。
  4. 如果计算区间已经足够小，则求解结束，得到最小值点，否则回到第二步。

对分法

对分法主要依赖 \(f'(x)\)，和求根的对分法基本类似。

算法：

  1. 起始的计算区间 \([x_l,x_r]=[a,b]\)。
  2. 选取中点 \(x_m = \frac12(x_l+x_r)\) ，计算 \(f'(x_m)\)。
  3. 如果中点 \(|f'(x_m)| \le \varepsilon\)，或者区间长度 \(|x_r-x_l|\le \varepsilon\) 已经足够小，则认为已经找到极值点，求解完成。
  4. 否则压缩区间：若 \(f'(x_m) > 0\)，取左侧为新的计算区间 \([x_l',x_r'] = [x_l,x_m]\) ；若 \(f'(x_m) < 0\)，取右侧为新的计算区间 \([x_l',x_r'] = [x_m,x_r]\)。

牛顿法

牛顿法需要依赖 \(f'(x),f''(x)\)，并且假设 \(f''>0\) 在计算区域内恒成立，否则难以保证收敛性。

本质是通过二阶泰勒展开将 \(f\) 在 \(x^{(k)}\) 局部近似为二次函数 \(q(h)\) ，然后更新到二次函数 \(q(h)\) 的极小值点，得到 \(x^{(k+1)}\)。

\[ q(h)= f(x^{(k)}) + f'(x^{(k)}) h + \frac{f''(x^{(k)})}2 h^2 \approx f(x^{(k)}+h) \]

极小值点 \[ h^* = -\frac{f'(x^{(k)})}{f''(x^{(k)})} \] 因此迭代格式为 \[ x^{(k+1)} = x^{(k)} + h^* = x^{(k)}-\frac{f'(x^{(k)})}{f''(x^{(k)})} \]

Remark： 关于多维的非线性方程最小值问题：自变量 \(X\in R^n\) ，目标函数 \(f: R^n \to R\) 足够光滑，假设目标函数存在极小值，那么牛顿迭代法的公式为 \[ X^{(k+1)} = X^{(k)} - H(X^{(k)})^{-1} \nabla f(X^{(k)}) \] 其中 \(\nabla f\) 是梯度向量，\(H = \nabla^2 f\) 是 Hesse 方阵。 \[ \begin{aligned} \nabla f(X)_i &= \frac{\partial f}{\partial X_i}(X)\\ H(X)_{ij} &= \nabla f(X)_{ij} = \frac{\partial^2 f}{\partial x_i \partial x_j}(X) \end{aligned} \] 要求 Hesse 方阵始终可逆，计算才可以继续，事实上我们还希望 Hesse 方阵始终对称正定，此时二阶近似才存在最小值点，否则可能不收敛。为了确保计算可以继续，从牛顿法发展出了拟牛顿法：用一族对称正定矩阵 \(H_k\) 来近似实际的 Hesse 方阵 \(H(X^{(k)})\) 参与计算。

非精确一维搜索

非精确一维搜索不要求直接求出最小值点，而是给定一个下降准则，求出一个点使得函数有足够的下降即可。

常用的下降准则有

Armijo准则：取定 \(\rho \in (0,1)\) ，要求 \(t\) 使得 \[ f(t) \le f(0) + \rho f'(0) t \] Goldstein准则：取定 \(\rho \in (0,\frac12)\) ，要求 \(t\) 使得 \[ \begin{aligned} f(t) &\le f(0) + \rho f'(0)t\\ f(t) &\ge f(0) + (1-\rho)f'(0) t \end{aligned} \] 注意：真正的最小值点也可能不满足Goldstein准则！

Wolfe准则：取定 \(\rho \in (0,\frac12)\) 以及 \(\sigma \in (\rho,1)\) ，要求 \(t\) 使得 \[ \begin{aligned} f(t) &\le f(0) + \rho f'(0)t \\ f'(t) &\ge \sigma f'(0) \end{aligned} \] 第二个条件通常会被加强为双边条件，称为强Wolfe准则 \[ |f'(t)| \le \sigma |f'(0)| = -\sigma f'(0) \]

基于这些准则可以设计很多具体的非精确一维搜索算法，主要的过程都是不断使用二分法、插值法或者步长乘以比例系数的方法得到新的试探点，判定准则是否成立，成立即可结束搜索，不成立则缩减搜索区间，继续尝试。

以下各有一个例子。

满足 Armijo 准则的回退法

基于Armijo准则设计的常用算法是回退法： 初始选取一个步长 \(t_0 = a\)，选取一个指数衰减因子 \(\gamma \in (0,1)\)，检查是否满足准则。 如果不满足准则，选取更小的步长 \(t_{k+1} = \gamma t_{k} = \gamma^{k} a\) 继续尝试，直到最终使用足够小的步长\(t\)，可以满足准则。

\[ t = \min_{k \in N}\{\gamma^k a : f(\gamma^k a) \le f(0) + \rho f'(0) \gamma^k a \} \]

满足 Goldstein 准则的对分法

取定初始的搜索区间 \([0,a]\)，以及两个参数 \(\rho \in (0,\frac12),\delta >0\)，计算\(f(0),f'(0)\) 用于准则的判定。

  1. 从\(t=0\)出发，当前搜索区间是 \([x_l,x_r]=[0,a]\)，选取中点作为试探点 \(x =\frac12 (x_l+x_r)\)，计算 \(f(x)\)
  2. 判定：验证第一条准则： \(f(x) \le f(0) + \rho f'(0)x\)
    1. 若满足第一条准则，验证第二条准则： \(f(x) \ge f(0) + (1-\rho)f'(0) x\)
      1. 若满足第二条准则，则搜索完成，输出 \(t=x\)
      2. 若不满足第二条准则，搜索区间改为右侧子区间 \([x_l',x_r'] = [x,x_r]\) ，下一次的试探点选取：如果右端点 \(x_r' < a\) 取中点 \(x' = \frac12(x_l'+x_r')\)，否则取 \(x' = (1+\delta)x\) 。继续判定
    2. 若不满足第一条准则，搜索区间改为左侧子区间 \([x_l',x_r'] = [x_l,x]\)，下一次的试探点仍然取作中点 \(x' = \frac12(x_l'+x_r')\) 。继续判定

满足 Wolfe 准则的算法

取定初始的搜索区间 \([0,a]\)，以及两个参数 \(\rho \in (0,\frac12),\sigma \in (\rho,1)\)，计算\(f(0),f'(0)\)用于后续判定。

  • 从\(t=0\)出发，当前搜索区间是 \([x_l,x_r]=[0,a]\)，选取适当的试探点 \(x \in (x_l,x_r)\)，计算 \(f(x)\)
  • 判定：验证第一条准则 \(f(x) \le f(0) + \rho f'(0)x\)
    1. 若满足第一条准则，验证第二条准则 \(f'(x) \ge \sigma f'(0)\)
      1. 若满足第二条准则，则搜索完成，输出 \(t=x\)
      2. 若不满足第二条准则，搜索区间改为右侧子区间 \([x_l',x_r'] = [x,x_r]\) ，下一次的试探点 \(x'\) 取为一个两点二次插值多项式的极小值点，继续判定
    2. 若不满足第一条准则，搜索区间改为左侧子区间 \([x_l',x_r'] = [x_l,x]\)，下一次的试探点 \(x'\) 取为一个两点二次插值多项式的极小值点，继续判定

更具体的细节：

当调整到左侧子区间\([x_l',x_r'] = [x_l,x]\)，取两点二次插值多项式 \(p(s),s\in [x_l,x]\)， 基于如下信息（我们还没有计算 \(f'(x)\) 的值） \[ \begin{aligned} p(x_l) &= f(x_l)\\ p'(x_l) &= f'(x_l)\\ p(x) &= f(x) \end{aligned} \] 得到 \(p(s)\) 的表达式 \[ \begin{aligned} p(s) &= p(x_l)\frac{(s-x)(s+x-2x_l)}{-(x_l-x)^2} + p'(x_l)\frac{(s-x_l)(s-x)}{x_l-x} + p(x)\frac{(s-x_l)^2}{(x-x_l)^2}\\ \end{aligned} \] 取新的试探点为极小值点 \[ x' = x_l +\frac12 \frac{(x_l-x)^2 p'(x_l)}{(p(x_l)-p(x))-(x_l-x)p'(x_l)} \]

当调整到右侧子区间\([x_l',x_r'] = [x,x_r]\)，取两点二次插值多项式 \(p(s),s\in [0,1]\)， 基于如下信息（我们已经计算了 \(f'(x)\) 的值） \[ \begin{aligned} p(x) &= f(x)\\ p'(x) &= f'(x)\\ p'(x_l) &= f'(x_l) \end{aligned} \] 得到 \(p(s)\) 的表达式 \[ \begin{aligned} p(s) &= p'(x_l)\frac{(s-x)^2}{2(x_l-x)} + p(x) + p'(x)\frac{(s-x)(s-2x_l+x)}{2(x-x_l)}\\ \end{aligned} \] 取新的试探点为极小值点 \[ x' = x - \frac{(x_l-x)p'(x)}{p'(x_l)-p'(x)} \]