学习要点 1. 理解优化算法的收敛性和收敛阶等概念. 2. 掌握一维牛顿法和割线法. 3. 掌握区间分割法. 4. 掌握几种常用的线搜索算法及其收敛性分析. 5.1 算法的收敛性与收敛速度 绝大多数非线性优化问题是无法求出其解析解的, 只能通过数值逼近的方法求其近 似解. 优化问题数值算法的基本原理是从某个初始点 x0 2 Rn 出发, 按照某种计算步骤产生 一列点 fxk : k = 0, 1, 2,    g, 使得 xk 逼近该优化问题的最优点 (或局部最优点). 当然, 最优点或局部最优点事先是不知道的, 有时候也用满足某种最优性条件的点代替它, 例如 KKT 点. 用 X. 表示要寻找的最优点的集合, fxkg 表示由某种算法产生的点列, 如果 d(xk,X.) := inf x∈X. kxk .. xk2 ! 0, k ! 1, (5.1) 则称 fxkg 弱收敛于 X.; 如果存在 x. 2 X. 使得 lim k→∞kxk .. x.k2 = 0, (5.2) 则称 fxkg 强收敛于 X.. 显然, 如果 fxkg 强收敛于 X., 则 fxkg 必弱收敛于 X.. 接下来讨论收敛速度的问题. 设 fxkg 收敛于 x. 2 X., 对于 p . 0, 称 Qp := lim k→∞ kxk+1 .. x.k2 kxk .. x.kp 2 (5.3) 为 fxkg 的商收敛因子, 简称为Q 因子; 称 Rp :=8<: lim k→∞kxk .. x.k1/k 2 , p = 1 lim k→∞kxk .. x.k1/pk 2 , p > 1 (5.4) 为 fxkg 的根收敛因子, 简称为R 因子. 称 OQ := inffp . 1 : Qp = 1g (5.5) 为 fxkg 的商收敛阶, 简称为Q 收敛阶; 称 OR := inffp . 1 : Rp = 1g (5.6) 为 fxkg 的根收敛阶, 简称为R 收敛阶. 对于收敛点列 fxkg, 可以证明其 Q 收敛阶不会超过其 R 收敛阶[34] . 目前文献中用得 比较多的是 Q 收敛阶, 以后我们提到的收敛因子都是指 Q 收敛因子, 收敛阶都是指 Q 收 敛阶. 定义 5.1 设点列 fxkg 收敛于 x., 如果 Q1 = 0, 则称 fxkg 超线性收敛于 x.; 如果 0 < Q1 < 1, 则称 fxkg 线性收敛于 x.; 如果 Q1 . 1, 则称 fxkg 次线性收敛于 x.. 如果 fxkg 超线性收敛于 x., 则有 lim k→∞ kxk+1 .. x.k2 kxk .. x.k2 = 0, (5.7) 由于 kxk+1 .. x.k2 kxk .. x.k2 . kxk+1 .. x.k2 kxk .. xk+1k2 + kxk+1 .. x.k2 = kxk+1 .. x.k2 kxk .. xk+1k2 1 + kxk+1 .. x.k2 kxk .. xk+1k2 , (5.8) 同理可得 kxk+1 .. x.k2 kxk .. xk+1k2 . kxk+1 .. x.k2 kxk .. x.k2 1 + kxk+1 .. x.k2 kxk .. x.k2 , (5.9) 因此式 (5.7) 成立当且仅当 lim k→∞ kxk+1 .. x.k2 kxk .. xk+1k2 = 0. (5.10) 这说明当相邻两次迭代的结果相差很小时, xk+1 已经距离 x. 很近了, 这就为算法提供了一 个很好的停机判断准则. 153 定义 5.2 设点列 fxkg 收敛于 x., 如果 Q2 = 0, 则称 fxkg 超二次收敛于 x.; 如果 0 < Q2 < 1, 则称 fxkg 二次收敛于 x.; 如果 Q2 = 1, 则称 fxkg 次二次收敛于 x.. 类似地, 可以定义超三次收敛、三次收敛、次三次收敛等概念. 5.2 一维牛顿法与割线法 本节我们考虑一维优化问题 min x∈R f(x) (5.11) 的优化算法. 我们假设 f(x) 具有二阶连续导数. 如果 x. 是式 (5.11) 的极小点, 则必有 f′(x.) = 0, 因此可以通过求方程 f′(x) = 0 的根求得候选极小点. 当然, 求这个方程的根也 并非易事, 需要用逐次逼近的方法求其近似值. 牛顿法 (也称为牛顿切线法) 就是这样一种 方法. 下面介绍牛顿法的基本思想. 首先选一个初始点 x0 2 R, 如果 f′(x0) = 0, 则 x0 就 是我们要找的点. 如果 f′(x0) 6= 0, 则可以考虑对 x0 作适当的修正, 例如使 x0 + h 比 x0 更接近 f′ 的根 x.. 如何确定修正量 h 呢? 我们当然希望 f′(x0 + h) = 0, 但这个方程与 f′(x) = 0 是等价的, 无法求出其精确解, 于是考虑在 x0 附近用一阶 Taylor 多项式逼近 f′, 即 f′(x0 + h) = f′(x0) + f′′(x0)h + o(jhj)  f′(x0) + f′′(x0)h, (5.12) 通过求解近似方程 f′(x0) + f′′(x0)h = 0 确定修正量 h. 这是一个线性方程, 很容易求得其 唯一的实根为 h = .. f′(x0) f′′(x0) , (5.13) 令 x1 = x0 + h = x0 .. f′(x0) f′′(x0) , (5.14) 则在一定的条件下, x1 确实比 x0 更接近 x.. 如果 f′(x1) 6= 0, 则可以重复以上过程对 x1 进行修正: 在 x1 点附近用一阶 Taylor 多 项式 f′(x1) + f′′(x1)h 近似代替 f′(x1 + h), 通过解线性方程 f′(x1) + f′′(x1)h = 0 确定修 正量 h = ..f′(x1)/f′′(x1), 令 x2 = x1 + h = x1 .. f′(x1)/f′′(x1). (5.15) 如果 f′(x2) 6= 0, 则再重复以上步骤, 直至某次迭代后得到的点 xk 满足 f′(xk) = 0 或者其 绝对值足够小, 停机并输出 xk. 以上计算过程可用算法表示如下. 154 Algorithm 1 (一维牛顿法 I) Input: The initial point x0 ∈ R; The tolerance bound ε; Output: The approximation of the extreme point x x ← x0; while |f′(x)| . ε do compute h = .f′(x)/f′′(x); x ← x + h; end while Output x; 关于以上算法的收敛性, 有下列结果. 定理 5.1 设 f 二次连续可微, f′(x.) = 0, f′′(x.) 6= 0, 则当 x0 充分靠近 x. 时, 有 xk ! x., 且 lim k→∞ jxk+1 .. x.j jxk .. x.j = 0, (5.16) 即算法 Algorithm 1 是超线性收敛的. 证明 根据算法的迭代公式得 xk+1 .. x.=xk .. f′(xk)/f′′(xk) .. x. =xk .. x. .. f′(x.) + ∫xk x. f′′(x)dx/f′′(xk) =∫xk x. [f′′(xk) .. f′′(x)] dx/f′′(xk), (5.17) 其中, 最后一个等号用到了条件 f′(x.) = 0. 由于 f′′(x.) 6= 0, 不妨设 f′′(x.) > 0; 由于 f′′ 连续, 因此存在 δ > 0 使得当 x, y 2 [x. .. δ, x. + δ] 时恒有 f′′(x) > f′′(x.)/2, jf′′(x) .. f′′(y)j < f′′(x.)/4, 从而有 jf′′(x) .. f′′(y)j f′′(x) < 1 2 , (5.18) 因此当 xk 2 [x. .. δ, x. + δ] 时有 jxk+1 .. x.j = ∫xk x. f′′(xk) .. f′′(x) f′′(xk) dx . 1 2jxk .. x.j, (5.19) 由此立刻推出 xk ! x.. 既然 xk ! x., 当 k ! 1 时必有 155 f′′(xk) .. f′′(x) f′′(xk) ! 0, 8 x 2 [x., xk], (5.20) 而且是一致收敛的, 因此有 jxk+1 .. x.j = ∫xk x. f′′(xk) .. f′′(x) f′′(xk) dx = jxk .. x.j  o(1), (5.21) 即式 (5.16) 成立. 算法 Algorithm 1 求出的只是导数的零点, 既可能是局部极小值点, 也可能是局部极大 值点, 一个完善的算法应能够判断是极小值点还是极大值点, 如果落入极大值点, 还要能够 跳出来继续寻找极小值点. 将以上因素都考虑进去的算法就是全局牛顿算法, 我们不作深 入讨论, 感兴趣的读者可参见文献 [31] 的 2.1 节. 一维牛顿法的关键迭代公式为 xk+1 = xk .. f′(xk) f′′(xk) , (5.22) 需要用到二阶导数 f′′(xk). 如果用差商近似代替 f′′(xk), 即 f′′(xk)  f′(xk) .. f′(xk.1) xk .. xk.1 , (5.23) 便得到了迭代公式 xk+1 = xk .. f′(xk)(xk .. xk.1) f′(xk) .. f′(xk.1) , (5.24) 由此得到的算法便是割线法, 具体算法如 Algorithm 2 所示. Algorithm 2 (割线法) Input: The initial point x0, x1 ∈ R; The tolerance bound ε; Output: The approximation of the extreme point xII xI ← x0; xII ← x1; while |f′(xII )| . ε do compute h = .f′(xII )(xII . xI)/(f′(xII ) . f′(xI )); xI ← xII ; xII ← xII + h; end while Output xII ; 156 关于割线法的收敛性的证明以及收敛阶的估计不作深入讨论, 感兴趣的读者可参见文 献 [31] 的 2.2 节. 5.3 区间分割法 本节讨论求解单峰函数最小值问题的方法——区间分割法. 称 f 是区间 [a, b] 上的单 峰函数 (unimodal function), 如果存在 x. 2 [a, b] 使得 f 在 [a, x.] 上严格单调递减, 在 [x., b] 上严格单调递增. 从单峰函数的定义不难发现, [a, b] 上的单峰函数 f 有唯一的全局极小值点 x.. 如何定 位极小值点 x. 呢? 我们可以在区间 [a, b] 上取两点 x1, x2 使得 a < x1 < x2 < b, 如果 f(x1) < f(x2), 则可以断定 x. 2 [a, x2]; 如果 f(x1) > f(x2), 则可以断定 x. 2 [x1, b]; 如果 f(x1) = f(x2), 则可以断定 x. 2 [x1, x2]. 综上所述, 通过计算两点的函数值 f(x1) 和 f(x2), 我们可以把 x. 定位在一个更小的区间中, 这就是区间分割法的核心思想. 通过反复 的区间分割, 可以给出极小值点 x. 任意精度的逼近. 如何选择区间分点 x1 和 x2 呢? 最容易想到的就是选择 x1 和 x2 为区间 [a, b] 的三等 分点: x1 = a + 1 3 (b .. a), x2 = a + 2 3 (b .. a), (5.25) 对应的区间分割算法就是算法 Algorithm 3. Algorithm 3 (区间分割算法 I) Input: The initial interval [a, b]; The tolerance bound ε; Output: The approximation of the extreme point xm am ← a; bm ← b; while bm . am . ε do x1 ← am + (bm . am)/3; x2 ← am + 2(bm . am)/3; if f(x1) < f(x2) then bm ← x2; else if f(x1) > f(x2) then am ← x1; else am ← x1; 157 bm ← x2; end if end if end while Output xm = (am + bm)/2; 算法 Algorithm 3 输出结果的精度由最终的小区间 [am, bm] 的长度决定, 满足不等式 jxm .. x.j . bm .. am 2 . (5.26) 因此算法的收敛速度取决于每次区间分割后极小值点 x. 所在的估计区间与原来的估计区 间相比缩小了多少. 设分割前的区间为 [am.1, bm.1], 分割后 x. 所在的区间为 [am, bm], 则 [am, bm] 或者是 [am.1, x2], 或者是 [x1, bm.1], 或者是 [x1, x2], 因此分割后与分割前的区间 长度之比满足 bm .. am bm.1 .. am.1 . maxfx2 .. am.1, bm.1 .. x1, x2 .. x1g bm.1 .. am.1 = maxfx2 .. am.1, bm.1 .. x1g bm.1 .. am.1 . (5.27) 当 x1 和 x2 是区间 [am.1, bm.1] 的三等分点时, 有 bm .. am bm.1 .. am.1 . 2 3 . (5.28) 经过 m 次区间分割后得到的近似极小值点 xm 的误差估计为 jxm .. x.j . bm .. am 2 . 1 2 2 3m (b .. a). (5.29) 为了提高计算效率, 我们可以在分割点的选择及算法的设计上做一些改进. 我们选择 分点 a < x < y < b 使得 x .. a = b .. y, y .. a b .. a = x .. a y .. a , (5.30) 令 γ = (y .. a)/(b .. a), 则有 γ = 1 γ .. 1, (5.31) 整理后得到一元二次方程 γ2 + γ .. 1 = 0, (5.32) 这个方程有两个实根, 但由于 γ > 0, 因此只能是 γ = p5 .. 1 2  0.618, (5.33) 158 这正是黄金分割比. 不难看出 b .. y b .. x = x .. a y .. a = γ. (5.34) 如果 f(x) < f(y), 则极小值点 x. 必落在区间 [a, y] 中, 只需令 x1 = a+(y..x), 则 a, x1, x, y 满足黄金分割条件, 只需再计算 f(x1) 的值便可进一步分割区间 [a, y], 更精确地定位 x.; 如果 f(x) . f(y), 则极小值点 x. 必落在区间 [x, b] 中, 只需令 y1 = b..(y..x), 则 x, y, y1, b 满足黄金分割条件, 只需再计算 f(y1) 的值便可进一步分割区间 [x, b], 更精确地定位 x.. 现在设计黄金分割算法如 Algorithm 4 所示. Algorithm 4 (黄金分割算法) Input: The initial interval [a, b]; The tolerance bound ε; Output: The approximation of the extreme point xm am ← a; bm ← b; γ ← (√5 . 1)/2; y ← am + γ(bm . am); x ← am + (bm . y); vx ← f(x); vy ← f(y); while bm . am . ε do if vx < vy then bm ← y; y ← x; vy ← vx; x ← am + (bm . y); vx ← f(x); else am ← x; x ← y; vx ← vy; y ← bm . (x . am); vy ← f(y); end if end while Output xm = (am + bm)/2; 黄金分割算法 Algorithm 4 每次分割后所得的估计区间的长度是分割前估计区间的长 度的 γ 倍, 因此有误差估计公式 jxm .. x.j . bm .. am 2 . 1 2 γm(b .. a), (5.35) 159 其中, m 是区间分割的次数. 与算法 Algorithm 3 相比有两点好处: 一是黄金分割比 γ < 2/3, 从而黄金分割算法收敛得更快; 二是黄金分割算法的每次迭代只需计算一个点的函数 值, 而算法 Algorithm 3 的每次迭代需要计算两个点的函数值, 因此黄金分割算法的效率 更高. 如果 f 存在一阶导数, 则极小值点 x. 是唯一满足 f′(x) = 0 的点, 因此可以通过求方 程 f′(x) = 0 的根来寻找 x.. 可以用对分法求方程 f′(x) = 0 的近似根: 设 c 是区间 [a, b] 的中点, 如果 f′(c) < 0, 则方程的根 x. 2 [c, b], 否则 x. 2 [a, c], 继续对分 [c, b] 或 [a, c] 可 得到更小的估计区间, 重复此过程可得到 x. 的任意精度的逼近. 具体算法如 Algorithm 5 所示. Algorithm 5 (对分法) Input: The initial interval [a, b]; The tolerance bound ε; Output: The approximation of the extreme point xm B ← 1; while b . a . ε do c ← (a + b)/2; if f′(c) < 0 then a ← c; else if f′(c) = 0 then Output xm = c; B ← 0; Stop; else b ← c; end if end if end while if B = 1 then Output xm = (a + b)/2; end if 对分法每次区间分割后估计区间的长度缩减为分割前估计区间的长度的一半, 因此有 误差估计公式 jxm .. x.j . 1 2m+1 (b .. a), (5.36) 其中, m 是区间分割的次数. 对分法的效率比三等分法和黄金分割法的都高, 前提是 f 可导. 160 5.4 线 搜 索 线搜索 (line search) 是指找多元函数在一条直线上的最大值点或最小值点. 线搜索 问题一般形式为 min s>0 f(xk + sdk), (5.37) 其中, xk 2 Rn 是某些优化算法迭代过程中到达的某一点, dk 2 Rn 代表接下来的行走方向, 通常满足 dTk rf(xk) < 0, 沿着这个方向走能够保证目标函数 f 的值下降. 线搜索问题之所 以重要是因为它出现在许多非线性优化算法的迭代步骤中, 对这些优算法的计算效率都有 直接影响. 由于线搜索问题本质上是一维优化问题, 因此可以用前面两节介绍的一维优化算法来 求解. 记 φ(s) = f(xk + sdk), 并假设 f 具有一阶连续偏导数, 则有 φ′(s) = dTk rf(xk + sdk), (5.38) 因此优化问题 (5.37) 在 s = s. 处取极小值的必要条件为 dTk rf(xk + s.dk) = 0. (5.39) 如果 f 是二阶连续可微的, 则有 φ′′(s) = dTk r2f(xk + sdk)dk, (5.40) 因此优化问题 (5.37) 在 s = s. 处取局部极小值的充分条件为 dTk rf(xk + s.dk) = 0, dTk r2f(xk + s.dk)dk > 0. (5.41) 如果 f 是凸函数, 则 φ 也是凸函数, 此时 φ′(s.) = 0 足以保证 φ 在 s = s. 处取得全局极 小值. 精确求解优化问题 (5.37) 称为精确线搜索 (exact line search), 前面两节介绍的牛 顿法、割线法和区间分割法都可以用于精确线搜索. 精确线搜索通常需要很多的计算时 间, 会严重拖慢优化算法的计算速度, 因此并不常用. 实践中用得更多的是非精确线搜索 (inexact line search), 目的是快速地找到某个近似解 sI 以满足某些指定的条件, 这些条 件能够保证非线性优化算法的整体收敛性. 一种比较常用的非精确线搜索是 Wolfe 线搜索, 它寻找 sI > 0 满足 f(x) .. f(x + sId) . ..sIb1dTrf(x), (5.42) dTrf(x + sId) . b2dTrf(x), (5.43) 其中, b1 和 b2 是两个常数, 且满足 0 < b1 . b2 < 1[35]. 条件 (5.42) 是为了保证当自变量 由 x 走到新的位置 x + sId 时目标函数 f 的值下降得足够多 (与步长 sI 成正比), 通常称 161