乘子罚函数
考虑等式约束问题:
min f ( x ) s.t. c i ( x ) = 0 , i = 1 , 2 , ⋯ , m \begin{aligned} \min ~~& f({x})\\ \text{s.t.} ~~& c_i(x) = 0, i = 1,2,\cdots,m \end{aligned} min s.t. f(x)ci(x)=0,i=1,2,⋯,m
其 Courant 罚函数
ϕ ( x , σ ) = f ( x ) + 1 2 σ ∑ i ( c i ( x ) ) 2 = f ( x ) + 1 2 σ c ( x ) T c ( x ) \begin{aligned} \phi(x,\sigma) = & f(x) + \frac{1}{2}\sigma\sum_i(c_i(x))^2 \\ = & f(x) + \frac{1}{2}\sigma c(x)^Tc(x) \end{aligned} ϕ(x,σ)==f(x)+21σi∑(ci(x))2f(x)+21σc(x)Tc(x)
在 σ k → ∞ \sigma_k \to \infty σk→∞ 时,得到一个局部极小点 x ∗ x^* x∗。确切地说,罚函数在极小点 x ( k ) x^{(k)} x(k) 处不再精确满足 c i ( x ) = 0 , i ∈ E c_i(x) = 0, i \in \mathcal{E} ci(x)=0,i∈E,而是被扰动为
c i ( k ) ≈ λ i ∗ σ k , i ∈ E c_i^{(k)} \approx \frac{
{\lambda_i}^*}{\sigma_k}, i \in \mathcal{E} ci(k)≈σkλi∗,i∈E
则有
c i ( k ) → 0 ∃ i , λ i ∗ ≠ 0 } ⇒ σ k → ∞ \left. \begin{aligned} c_i^{(k)} \to 0 \\ \exists i,~~ \lambda_i^* \neq 0 \end{aligned} \right \} \Rightarrow \sigma_k \to \infty ci(k)→0∃i, λi∗=0}⇒σk→∞
动机
构造新的罚函数,使得固定某罚参数后,无约束优化问题的解与原始问题的相同!
改造约束函数
对 c i ( x ) c_i(x) ci(x) 进行平移,即在罚函数中使用 c i ( x ) − θ i c_i(x) - \theta_i ci(x)−θi 代替 c i ( x ) c_i(x) ci(x),使得对于有限的 σ \sigma σ, ϕ \phi ϕ 可以在 x ∗ x^* x∗ 处取得极小值。该种方法由 Powell 于 1969 年提出。
ϕ ( x , θ , σ ) = f ( x ) + 1 2 σ ∑ i ∈ E ( c i ( x ) − θ i ) 2 \phi(x,\theta,\sigma) = f(x) + \frac{1}{2}\sigma \sum_{i \in \mathcal{E}} (c_i(x) - \theta_i)^2 ϕ(x,θ,σ)=f(x)+21σi∈E∑(ci(x)−θi)2
例 考虑问题
KaTeX parse error: No such environment: equation* at position 8: \begin{̲e̲q̲u̲a̲t̲i̲o̲n̲*̲}̲\begin{aligned}…
的乘子罚函数
ϕ ( x , σ ) = x + 1 2 σ ( 1 − x − θ ) 2 \phi(x,\sigma) = x + \frac{1}{2}\sigma(1 - x - \theta)^2 ϕ(x,σ)=x+21σ(1−x−θ)2
如下图所示:
由上图可知,该方法能够保证在 σ \sigma σ 有限的情况下,通过调整 θ \theta θ 的大小以保证 x ∗ x^* x∗ 就是 ϕ ( x , θ , σ ) \phi(x,\theta,\sigma) ϕ(x,θ,σ) 的最优解,从而避免了极限 σ → ∞ \sigma \to \infty σ→∞ 的病态问题。
改造目标函数
改造目标函数以避免对系统的扰动,即对适中的 σ k \sigma_k σk,近似极小点更好地满足等式约束 c i ( x ) = 0 c_i(x) = 0 ci(x)=0。该方法由 Hestenes 于 1969 年提出,通过在罚函数中引入 Lagrange 乘子的显式估计可达到这种目标。
ϕ ( x , λ , σ ) = f ( x ) + λ T c ( x ) + 1 2 σ c ( x ) T c ( x ) \phi(x,\lambda,\sigma) = f(x) + \lambda^T c(x) + \frac{1}{2} \sigma c(x)^Tc(x) ϕ(x,λ,σ)=f(x)+λTc(x)+21σc(x)Tc(x)
该方法由于引入了乘子项 λ T c ( x ) \lambda^Tc(x) λTc(x),因此也称为乘子罚函数;又因为在 Lagrange 函数中引入罚项 1 2 σ c ( x ) T c ( x ) \frac{1}{2} \sigma c(x)^Tc(x) 21σc(x)Tc(x) 也可以得到,因而又称增广 Lagrange 函数。
Powell-Hestenes 函数
其实,令
λ i = − θ i σ , i = 1 , 2 , ⋯ , m \lambda_i = -\theta_i\sigma, i = 1,2,\cdots,m λi=−θiσ,i=1,2,⋯,m
Powell 函数展开后只比 Hestenes 函数多了一个与 x x x 无关项 1 2 σ ∑ i θ i 2 \frac{1}{2}\sigma\sum_i\theta_i^2 21σ∑iθi2,因此将 Hestenes 函数也称为 Powell-Hestenes 函数。
乘子罚函数的性质及特点
精确性
定理 1 如果 x ∗ x^* x∗, λ ∗ \lambda^* λ∗ 处二阶充分条件成立,则 ∃ σ ′ ≥ 0 \exists \sigma'\geq 0 ∃σ′≥0,对于 ∀ σ ≥ σ ′ \forall \sigma \geq \sigma' ∀σ≥σ′, x ∗ x^* x∗ 为 ϕ ( x , λ ∗ , σ ) \phi(x,\lambda^*,\sigma) ϕ(x,λ∗,σ) 的严格局部极小点,即 x ∗ = x ( λ ∗ ) x^* = x(\lambda^*) x∗=x(λ∗)。
在实际中,并不能得到 λ ∗ \lambda^* λ∗ 的精确值,该结论表明:若 λ \lambda λ 是 λ ∗ \lambda^* λ∗ 的好的估计,那么对于不是很大的 σ \sigma σ,通过极小化 ϕ ( x , λ , σ ) \phi(x,\lambda,\sigma) ϕ(x,λ,σ) 可以得到 λ ∗ \lambda^* λ∗ 的好的估计。
乘子的特点
- 子问题是光滑的
- 可以利用"使用导数的方法"求解子问题
- 一定条件下,不需要 σ k → ∞ \sigma_k \to \infty σk→∞
- 避免了病态的 Hessian 矩阵
算法
选取充分大的控制参数 σ \sigma σ 后,用 λ \lambda λ 作为序列极小化算法中的控制参数,有如下算法:
算法 1 乘子法求解算法
-
选择一个乘子序列 { λ ( k ) } \{\lambda^{(k)}\} { λ(k)},使得 λ ( k ) → λ ∗ \lambda^{(k)} \to \lambda^* λ(k)→λ∗
-
do
对于 ∀ λ ( i ) \forall \lambda^{(i)} ∀λ(i),求解
min x ϕ ( x , λ ( k ) , σ ) \min_x~~\phi(x,\lambda^{(k)},\sigma) xmin ϕ(x,λ(k),σ)
的最优解 x ( λ ( k ) ) x(\lambda^{(k)}) x(λ(k))。while c ( x ( λ ( k ) ) ) c(x(\lambda^{(k)})) c(x(λ(k))) 充分小
实际中,如何确定乘子序列,以及罚参考如何选择将是需要解决的首要问题。
乘子更新
欲求解 ϕ ( x , λ , σ ) \phi(x,\lambda,\sigma) ϕ(x,λ,σ) 的极小点
x ( λ ) = arg min x ∈ R n ϕ ( x , λ , σ ) x(\lambda) = \mathop{\arg\min}_{x \in \mathbb{R}^n} \phi(x,\lambda,\sigma) x(λ)=argminx∈Rnϕ(x,λ,σ)
需求
∇ ϕ ( x , λ , σ ) = ∇ f ( x ) + ∑ i ∈ E ( λ i + σ c i ( x ) ) ∇ c i ( x ) = 0 \begin{aligned} \nabla\phi(x,\lambda,\sigma) &= \nabla f(x) + \sum_{i \in \mathcal{E}}(\lambda_i + \sigma c_i(x)) \nabla c_i(x) \\ &= 0 \end{aligned} ∇ϕ(x,λ,σ)=∇f(x)+i∈E∑(λi+σci(x))∇ci(x)=0
在确定 x ( λ ) x(\lambda) x(λ) 后,得
ψ ( λ ) : = ϕ ( x ( λ ) , λ , σ ) \psi(\lambda) \mathop{\colon =} \phi(x(\lambda),\lambda,\sigma) ψ(λ):=ϕ(x(λ),λ,σ)
则
λ ∗ = arg max λ ∈ R m ψ ( λ ) \lambda^* = \mathop{\arg\max}_{\lambda\in \mathbb{R}^m}\psi(\lambda) λ∗=argmaxλ∈Rmψ(λ)
计算可得
∇ ψ ( λ ) = c ( x ( λ ) ) ∇ 2 ψ ( λ ) = d c d λ = − A T W σ − 1 A ∣ x ( λ ) \begin{aligned} \nabla \psi(\lambda) &=& c(x(\lambda)) \\ \nabla^2 \psi(\lambda) &=& \frac{d c}{d \lambda} = - A^TW_{\sigma}^{-1}A |_{x(\lambda)} \end{aligned} ∇ψ(λ)∇2ψ(λ)==c(x(λ))dλdc=−ATWσ−1A∣x(λ)
其中
W σ = ∇ x 2 ϕ ( x , λ , σ ) = W + σ A A T W = ∇ 2 f ( x ) + ∑ i [ λ i + σ c i ( x ) ] ∇ 2 c i ( x ) \begin{aligned} W_{\sigma} &=& \nabla^2_x\phi(x,\lambda,\sigma) = W + \sigma AA^T \\ W &=& \nabla^2f(x) + \sum_{i} [\lambda_i + \sigma c_i(x)] \nabla^2c_i(x) \end{aligned} WσW==∇x2ϕ(x,λ,σ)=W+σAAT∇2f(x)+i∑[λi+σci(x)]∇2ci(x)
且 A = [ ∇ c 1 ( x ) , c 2 ( x ) , ⋯ , c m ( x ) ] A = [\nabla c_1(x), c_2(x),\cdots,c_m(x)] A=[∇c1(x),c2(x),⋯,cm(x)]。
牛顿法
按照牛顿法迭代选取序列 { λ ( k ) } \{\lambda^{(k)}\} {
λ(k)},确定初始 λ 0 \lambda^0 λ0,再利用
λ ( k + 1 ) = λ ( k ) + ( A T W σ − 1 A ) − 1 c ∣ x ( λ ( k ) ) \lambda^{(k +1)} = \lambda^{(k)} + (A^T W_{\sigma}^{-1}A)^{-1}c|_{x(\lambda^{(k)})} λ(k+1)=λ(k)+(ATWσ−1A)−1c∣x(λ(k))
这种方式用到了 W σ W_{\sigma} Wσ 二阶导数的显式表达,只有一阶导可用时,可以利用拟牛顿法进行近似。
对于充分大的 σ \sigma σ,有
( A T W σ − 1 A ) − 1 ≈ σ I (A^T W_{\sigma}^{-1}A)^{-1} \approx \sigma I (ATWσ−1A)−1≈σI
因此,得到迭代
λ ( k + 1 ) = λ ( k ) + σ c i ( k ) , i ∈ E \lambda^{(k +1)} = \lambda^{(k)} + \sigma c_i^{(k)}, i \in \mathcal{E} λ(k+1)=λ(k)+σci(k),i∈E
该式不需要任何导数,通过充分大的 σ \sigma σ 可以让 λ ( k ) \lambda^{(k)} λ(k) 以任意快的速度线性收敛于 λ ∗ \lambda^* λ∗。
例 考虑等式规划问题
KaTeX parse error: No such environment: equation* at position 8: \begin{̲e̲q̲u̲a̲t̲i̲o̲n̲*̲}̲\begin{aligned}…
结果如下表所示:
罚参数更新
给定初始罚参数 σ \sigma σ,解得最优解 x ( 0 ) x^{(0)} x(0) 和 乘子的估计 λ ( 1 ) \lambda^{(1)} λ(1)
第 k k k 次迭代,固定参数 σ \sigma σ 和 λ ( k ) \lambda^{(k)} λ(k),得到 Warm-start 技术
x ′ = arg min x ∈ R n ϕ ( x , λ ( k ) , σ ) x' = \mathop{\arg\min}_{x\in \mathbb{R}^n} \phi(x,\lambda^{(k)},\sigma) x′=argminx∈Rnϕ(x,λ(k),σ)
根据需要更新的罚参数
-
若 ∥ c ( x ′ ) ∥ ∞ > 1 4 ∥ c ( k − 1 ) ∥ ∞ \|c(x')\|_{\infty} > \frac{1}{4} \|c^{(k-1)}\|_{\infty} ∥c(x′)∥∞>41∥c(k−1)∥∞,令
σ = 10 σ , x ( k − 1 ) = x ′ \sigma = 10 \sigma,x^{(k-1)} = x' σ=10σ,x(k−1)=x′
重复第 k k k 次迭代 -
否则,令
x ( k ) = x ′ , λ ( k + 1 ) = λ ( k ) + σ c ( x ( k ) ) x^{(k)} = x', \lambda^{(k+1)} = \lambda^{(k)} + \sigma c(x^{(k)}) x(k)=x′,λ(k+1)=λ(k)+σc(x(k)) -
令 k = k + 1 k = k +1 k=k+1,进行第 k k k 次迭代。
参考资料
[1] 刘红英,夏勇,周永生. 数学规划基础,北京,北京航空航天大学出版社,2012.