Numerical Optimization Ch17. Penalty and Augmented Lagragian Methods

第十七章: 惩罚函数法与增广Lagrange函数法

一些重要的求解约束优化的方法 将原问题替换为一系列子问题, 在子问题中原本的约束被替换成加在目标函数上的项. 本章我们介绍属于此类的三种方法.

• 二次罚函数法the quadratic penalty method是在目标函数上增加每个约束违反度平方的某个倍数. 这种方法比较简单、直观. 尽管它有许多的缺陷, 但实际中还是被经常使用.
• 非光滑精确罚函数法the nonsmooth exact penalty methods是用一个(而不是一系列)无约束问题替代原本的约束问题. 使用这种罚函数, 我们可以仅利用一次无约束极小化找到解, 但随之而来, 非光滑性可能会制造一些麻烦. 此类的常见罚函数比如$l_1$罚函数.
• 乘子法the methods of multipliers或增广Lagrange函数法the augmented Lagrangian method是另一种精确罚函数法, 其中显式估计了Lagrange乘子, 避免了二次罚函数带来的病态问题.
• 对数障碍函数法log-barrier method使用对数项来防止迭代点过于接近可行域边界. 此法是非线性规划内点法的部分基础. 我们放在第十九章讲述.

1. 二次罚函数法

1.1 动机

考虑使用单一个函数代替约束优化问题, 其中那个函数由以下组成:

• 约束优化问题的目标函数, 加上
• 对每个约束的一个附加项, 其在当前点违反约束时为正, 否则取0.

许多方法都会通过定义一系列这样的罚函数求解问题, 其中罚项乘上了一个正系数(惩罚因子). 此系数越大, 我们对违反约束的现象越不能容忍, 从而越强制罚函数的极小点靠近约束问题的可行域.

此类最简单的罚函数为二次罚函数(又称Courant罚函数), 其中的罚项记为约束违反度的平方. 我们先在等式约束问题下讨论这一方法.$$\underset{x}{\min }f(x),\mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}{c}_i(x)=0,i\in \mathcal{E}.$$对此, 二次罚函数为$$Q(x;\mu )\stackrel{def}{=}f(x)+\frac{\mu }{2}\sum _{i\in \mathcal{E}}{c}_i^2(x),$$其中$\mu >0$为惩罚因子(penalty parameter). 当$\mu$趋向于$\infty$, 我们就对约束违反惩罚得愈加猛烈. 直观上, 我们很容易想到要构作一列 { μ k } : μ ↑ ∞ , k → ∞ \{\mu_k\}:\mu\uparrow\infty,k\to\infty {μk​}:μ↑∞,k→∞, 并对每个$k$求得$Q(x;{\mu }_k)$的近似极小点$x_k$. 由于上述罚函数是光滑的, 因此我们可以使用无约束优化的技术求$x_k$. 在搜索$x_k$时, 我们可以使用先前的极小点$x_{k-1},x_{k-2}$作为算法初始的迭代点. 对于合适的 { μ k } \{\mu_k\} {μk​}和初始点, 我们在每次无约束极小过程中都无需消耗过多的计算量.

例1 考虑以下等式约束问题$$\min x_1+x_2,\mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}{x}_1^2+x_2^2-2=0,$$其解为$(-1,-1{)}^T$且二次罚函数为$$Q(x;\mu )=x_1+x_2+\frac{\mu }{2}(x_1^2+x_2^2-2{)}^2.$$下图为$\mu =1$时, $Q$的等高线图. 从图中可见$Q$有一差不多是$(-1.1,-1.1{)}^T$的极小点, 也有一接近$(0.3,0.3{)}^T$的极大点.

而下图则对应$\mu =10$, 因此不再可行域$x_1^2+x_2^2=2$的点会遭受比上图更大的惩罚. 此图中$Q$的极小点愈加靠近解. 同时还有一局部极大点在$(0,0{)}^T$附近, 且$Q$在出可行域后就迅速地趋向于$\infty$.

但实际情形一般没有例1中那么良态. 对于一给定的$\mu$, 即使原约束问题有唯一解, 惩罚函数也可能是下无界的. 例如,$$\min -5x_1^2+x_2^2,\mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}{x}_1=1,$$其解为$(1,0{)}^T$. 而对此任何小于10的惩罚因子, 其罚函数都下无界. 不过, 这点对于本章讨论的所有惩罚函数都是存在的.

对于一般的约束优化问题$$\underset{x}{\min }f(x),\mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}{c}_i(x)=0,i\in \mathcal{E},c_i(x)\ge 0,i\in \mathcal{I}.$$定义二次罚函数为$$Q(x;\mu )\stackrel{def}{=}f(x)+\frac{\mu }{2}\sum _{i\in \mathcal{E}}{c}_i^2(x)+\frac{\mu }{2}\sum _{i\in \mathcal{I}}{\left([c_i(x){]}^{-}\right)}^2,$$这里$[y{]}^{-}=\max (-y,0)$. 此时$Q$没有光只有等式约束的惩罚函数光滑, 也不如目标函数和约束函数光滑. 例如若有一个不等式约束为$x_1\ge 0$, 则$\min (0,x_1{)}^2$就没有连续的二阶导数, $Q$就不再二阶连续可微.

1.2 算法框架

基于二次罚函数的一般算法框架如下.

框架1 (Quadratic Penalty Method)
Given ${\mu }_0>0$, a nonnegative sequence { τ k } \{\tau_k\} {τk​} with ${\tau }_k\to 0$, and a starting point $x_0^s$;
for $k=0,1,2,\dots$
\quad\quad Find an approximate minimizer $x_k$ of $Q(\cdot ;{\mu }_k)$, starting at $x_k^s$, and terminating when $\Vert {\nabla }_{x}Q(x;{\mu }_k)\Vert \le {\tau }_k$;
\quad\quad if final convergence test satisfied
\quad\quad\quad\quad stop with approximate solution $x_k$;
\quad\quad end (if)
\quad\quad Choose new penalty parameter ${\mu }_{k+1}>{\mu }_k$;
\quad\quad Choose new starting point $x_{k+1}^s$;
end (for)

我们可以基于在每次迭代极小化罚函数的困难度自适应地选取惩罚因子序列 { μ k } \{\mu_k\} {μk​}:

• 当极小化$Q(x;{\mu }_k)$对某一$k$来说比较昂贵, 我们就选取比${\mu }_k$稍微大一点的${\mu }_{k+1}$, 比如${\mu }_{k+1}=1.5{\mu }_k$.
• 若我们比较容易就能找到$Q(x;{\mu }_k)$的近似极小点, 我们就尝试更大幅度的增长, 比如${\mu }_{k+1}=10{\mu }_k$.

框架1的收敛理论给予了非负序列 { τ k } \{\tau_k\} {τk​}选取充分的自由度: 仅需${\tau }_k\to 0$, 以保证算法随着迭代, 其精度在不断提升.

正如之前讨论过的, 我们无法保证$\Vert {\nabla }_{x}Q(x;{\mu }_k)\Vert \le {\tau }_k$一定会成立. 因此实际的实施必须包括当约束违反度下降不够快, 或者当迭代点趋向于发散时, 增大惩罚因子(或存储初始点)的防护措施.

当只有等式约束时, $Q(x;{\mu }_k)$光滑, 因此可使用无约束极小化的算法求得近似解$x_k$. 不过, 当${\mu }_k$变大, $Q(x;{\mu }_k)$的极小化会变得愈发困难(求解系统的条件数变大), 除非我们使用特殊的手段计算搜索方向.

一方面, Hessian阵${\nabla }_{xx}^{2}Q(x;{\mu }_k)$会在极小点附近变得愈发病态. 这一点就足以使得许多无约束算法(如拟牛顿法、共轭梯度法)效果变差. 而诸如牛顿法这般对Hessian条件数不敏感的算法, 则可能会遭遇较大${\mu }_k$带来的另外两个问题:

1. 在我们求解线性系统计算牛顿步时, 病态的${\nabla }_{xx}^{2}Q(x;{\mu }_k)$可能会引发数值不稳定. 后面我们会指出, 这种影响并不严重, 我们可以重构牛顿方程.
2. 即使$x$很接近$Q(\cdot ;{\mu }_k)$的极小点, 但对于$Q(x;{\mu }_k)$的二阶Taylor展开仅在$x$的一个小邻域内才是合适的. 这一点可从第二张图看出来, 其中$Q$在极小点附近的等高线呈现出一种"香蕉"状而不是二次函数典型的椭圆型. 由于牛顿法是基于二次模型的, 因此其产生的迭代步可能无法快速收敛到$Q(x;{\mu }_k)$的极小点. 对此, 我们需要谨慎地选择初始点$x_{k+1}^s$或直接设置$x_{k+1}^s=x_k$, 并选取${\mu }_{k+1}$只比${\mu }_k$大一点.

1.3 二次罚函数法的收敛性

我们在下面两个定理里介绍二次罚函数方法的一些收敛性质. 我们将讨论限制在等式约束问题.

对于第一个结果, 我们假设对每个${\mu }_k$, 罚函数$Q(x;{\mu }_k)$都有(有限个)极小点.

定理1 设$x_k$为框架1中$Q(x;{\mu }_k)$的精确全局极小点, 且${\mu }_k\uparrow \infty$. 则序列 { x k } \{x_k\} {xk​}的每个聚点都是原问题的全局解.

证明: 令$\bar{x}$为原问题的全局解. 即有$$f(\bar{x})\le f(x),\forall x:c_i(x)=0,i\in \mathcal{E}.$$由于对每个$k$, $x_k$极小化$Q(\cdot ;{\mu }_k)$, 于是$Q(x_k;{\mu }_k)\le Q(\bar{x};{\mu }_k)$, 这就得到不等式$$f(x_k)+\frac{{\mu }_k}{2}\sum _{i\in \mathcal{E}}{c}_i^2(x_k)\le f(\bar{x})+\frac{{\mu }_k}{2}\sum _{i\in \mathcal{E}}{c}_i^2(\bar{x})=f(\bar{x}).$$重新整理可得$$\sum _{i\in \mathcal{E}}{c}_i^2(x_k)\le \frac{2}{{\mu }_k}[f(\bar{x})-f(x_k)].$$假设$x^{\ast }$为 { x k } \{x_k\} {xk​}的一个聚点, 因此存在 { 1 , 2 , … , n , … } \{1,2,\ldots,n,\ldots\} {1,2,…,n,…}的无穷子列$\mathcal{K}$使得$$\underset{k\to \mathcal{K}}{\lim }{x}_k=x^{\ast }.$$对上面不等式两边取极限$k\to \infty ,k\in \mathcal{K}$, 我们有$$\sum _{i\in \mathcal{E}}{c}_i^2(x^{\ast })=\underset{k\in \mathcal{K}}{\lim }\sum _{i\in \mathcal{E}}{c}_i^2(x_k)\le \underset{k\in \mathcal{K}}{\lim }\frac{2}{{\mu }_k}[f(\bar{x})-f(x_k)]=0.$$因此$c_i(x^{\ast })=0,i\in \mathcal{E}$, 从而$x^{\ast }$可行. 进一步, $$f(x^{\ast })\le f(x^{\ast })+\underset{k\in \mathcal{K}}{\lim }\frac{{\mu }_k}{2}\sum _{i\in \mathcal{E}}{c}_i^2(x_k)\le f(\bar{x}).$$由于$x^{\ast }$可行, 且其函数值不大于全局解$\bar{x}$处的函数值, 因此$x^{\ast }$也是全局解. 证毕.

此结论对带不等式的问题也成立. 但由于需要我们求每个子问题的全局极小点, 因此此性质一般并不能成立. 下一结论则允许我们不精确(但精确度要提升)极小化$Q(\cdot ;{\mu }_k)$. 相较于定理1, 它表示 { x k } \{x_k\} {xk​}可能会收敛到不可行点或是KKT点. 它也说明了在一定情形下, $-{\mu }_k{c}_i(x_k)$可以当做Lagrange乘子${\lambda }_i^{\ast }$的估计. 这一点对于我们在第3节讨论增广Lagrange函数法时比较重要.

为建立定理, 我们需乐观地假设对所有$k$, $\Vert {\nabla }_{x}Q(x;{\mu }_k)\Vert \le {\tau }_k$都会成立.

定理2 设框架1中的容忍限序列和惩罚因子序列满足${\tau }_k\to 0,{\mu }_k\uparrow \infty$, 则若序列 { x k } \{x_k\} {xk​}的聚点$x^{\ast }$不可行, 它就是函数$\Vert c(x){\Vert }^2$的稳定点. 若聚点$x^{\ast }$可行且约束梯度$\nabla c_i(x^{\ast })$线性无关, 则$x^{\ast }$为原问题的KKT点. 于此, 对任一收敛到$x^{\ast }$的子列, 即$\forall \mathcal{K}:{\lim }_{k\in \mathcal{K}}{x}_k=x^{\ast }$, 我们有$$\underset{k\in \mathcal{K}}{\lim }-{\mu }_k{c}_i(x_k)={\lambda }_i^{\ast },\forall i\in \mathcal{E},$$这里${\lambda }^{\ast }$为满足等式约束原问题KKT条件的乘子.

证明: 对$Q(x;{\mu }_k)$求导, 得到$${\nabla }_{x}Q(x_k;{\mu }_k)=\nabla f(x_k)+\sum _{i\in \mathcal{E}}{\mu }_k{c}_i(x_k)\nabla c_i(x_k),$$于是从终止条件, 我们有$$\Vert \nabla f(x_k)+\sum _{i\in \mathcal{E}}{\mu }_k{c}_i(x_k)\nabla c_i(x_k)\Vert \le {\tau }_k.$$利用三角不等式, 有$$\Vert \sum _{i\in \mathcal{E}}{c}_i(x_k)\nabla c_i(x_k)\Vert \le \frac{1}{{\mu }_k}[{\tau }_k+\Vert \nabla f(x_k)\Vert ].$$令$x^{\ast }$为迭代点序列的聚点. 于是存在无穷指标子列$\mathcal{K}$使得${\lim }_{k\in \mathcal{K}}{x}_k=x^{\ast }$. 对上式取$k\to \infty ,k\in \mathcal{K}$的极限, 可得$$\sum _{i\in \mathcal{E}}{c}_i(x^{\ast })\nabla c_i(x^{\ast })=0.$$我们可能有$c_i(x^{\ast })\ne 0$(若约束梯度$\nabla c_i(x^{\ast })$线性相关), 但这起码说明$x^{\ast }$为$\Vert c(x){\Vert }^2$的稳定点.

若$\nabla c_i(x^{\ast })$线性无关, 则$c_i(x^{\ast })=0,\forall i\in \mathcal{E}$, 因此$x^{\ast }$可行. 于是KKT条件之可行性条件成立. 记约束的Jacobi阵为$A$,向量$-{\mu }_{k}c(x_k)$为${\lambda }_k$, 则有$$A(x_k{)}^T{\lambda }_k=\nabla f(x_k)-{\nabla }_{x}Q(x_k;{\mu }_k),\Vert {\nabla }_{x}Q(x_k;{\mu }_k)\Vert \le {\tau }_k.$$对$k\in \mathcal{K}$充分大, $A(x_k)$行满秩, 于是$A(x_k)A(x_k{)}^T$非奇异. 上式左乘$A(x_k)$并整理可得$${\lambda }_k={\left[A(x_k)A(x_k{)}^T\right]}^{-1}A(x_k)[\nabla f(x_k)-{\nabla }_{x}Q(x_k;{\mu }_k)].$$取极限$k\to \infty ,k\in \mathcal{K}$, 可得$$\underset{k\in \mathcal{K}}{\lim }{\lambda }_k={\lambda }^{\ast }={\left[A(x^{\ast })A(x^{\ast }{)}^T\right]}^{-1}A(x^{\ast })\nabla f(x^{\ast }).$$又由前面取极限可知$$\nabla f(x^{\ast })-A(x^{\ast }{)}^T{\lambda }^{\ast }=0,$$于是${\lambda }^{\ast }$满足KKT条件之稳定性条件. 所以, $x^{\ast }$为原问题的KKT点, 且由唯一的Lagrange乘子${\lambda }^{\ast }$. 证毕.

需重复强调的是, 若聚点$x^{\ast }$不可行, 则它至少是$\Vert c(x){\Vert }^2$的稳定点. 牛顿类的算法总是会陷在这种类型的点. 这有点像我们在第十一章中讨论非线性方程组的情形. 在原问题不可行时, 我们往往可以观察到二次罚函数算法收敛于$\Vert c(x){\Vert }^2$的稳定点或极小点.

若考虑不等式约束, 则$x^{\ast }$为函数$\Vert [c(x){]}^{-}{\Vert }^2$的稳定点, 其中$[c(x){]}^{-}$定义为$$[c(x){]}_i^{-}=\{\begin{array}{ll}{c}_i(x),& i\in \mathcal{E},\\ [c_i(x){]}^{-},& i\in \mathcal{I}.\end{array}$$但所谓的KKT点不一定成立, 因为此处无法估计乘子的正负.

1.4 病态与重构

我们现在来验证Hessian阵${\nabla }_{xx}^{2}Q(x;{\mu }_k)$的病态内在. 对这一矩阵和在其他罚函数和障碍函数法中出现的Hessian阵性质的理解, 将有助于我们选取和设计合适、高效的算法.

这里的Hessian阵为$${\nabla }_{xx}^{2}Q(x;{\mu }_k)={\nabla }^{2}f(x)+\sum _{i\in \mathcal{E}}{\mu }_k{c}_i(x){\nabla }^2{c}_i(x)+{\mu }_{k}A(x{)}^{T}A(x).$$当$x$接近$Q(\cdot ;{\mu }_k)$的极小点且定理2的条件满足时, 我们从定理2的结论可知, 上式右端前两项就差不多是Lagrange函数的Hessian. 具体地, $${\nabla }_{xx}^{2}Q(x;{\mu }_k)\approx {\nabla }_{xx}^2\mathcal{L}(x,{\lambda }^{\ast })+{\mu }_{k}A(x{)}^{T}A(x).$$依此, ${\nabla }_{xx}^{2}Q(x;{\mu }_k)$就接近于以下两个部分的和:

• (Lagrange项)独立于${\mu }_k$的矩阵; 和
• 秩为$|\mathcal{E}|$, 且非零特征值与${\mu }_k$同阶的矩阵.

由于约束的数量$|\mathcal{E}|$通常小于$n$, 因此第二项奇异. 因此整个矩阵的某些特征值会趋近常数, 而其他的则与${\mu }_k$同阶. 因${\mu }_k\to \infty$, 通过条件数可知, ${\nabla }_{xx}^{2}Q(x;{\mu }_k)$会变得愈发病态.

Hessian病态的一个后果是, 以下牛顿步的计算可能会不够精确:$${\nabla }_{xx}^{2}Q(x;{\mu }_k)p=-{\nabla }_{x}Q(x;{\mu }_k).$$一般来说, 不论使用何种直接计算手段, 此病态系统都会导致在$p$上的巨大计算误差. 基于同样的原因, 除非有预处理的策略, 否则迭代计算的方法也不会有较好的表现.

不过, 我们可以重构以上牛顿方程以避免病态. 引入一个新的变量$\zeta :=\mu A(x)p$, 我们会发现满足牛顿方程的$p$同样满足以下系统:$$\left[\begin{array}{cc}{\nabla }^{2}f(x)+\sum _{i\in \mathcal{E}}{\mu }_k{c}_i(x){\nabla }^2{c}_i(x)& A(x{)}^T\\ A(x)& -(1/{\mu }_k)I\end{array}\right]\left[\begin{array}{c}p\\ \zeta \end{array}\right]=\left[\begin{array}{c}-{\nabla }_{x}Q(x;{\mu }_k)\\ 0\end{array}\right].$$当$x$离$x^{\ast }$较近时, 此系统的系数矩阵并不会有与${\mu }_k$同阶的大特征值. 事实上, 当${\mu }_k\to \infty$, 可将此系统的系数矩阵近似看做$$\left[\begin{array}{cc}{\nabla }_{xx}^2\mathcal{L}(x,{\lambda }^{\ast })& A(x{)}^T\\ A(x)& O\end{array}\right],$$进一步若假设${\nabla }_{xx}^2\mathcal{L}(x,{\lambda }^{\ast })$非奇异, 利用分块矩阵的初等变换, 可得$$\left[\begin{array}{cc}{\nabla }_{xx}^2\mathcal{L}(x,{\lambda }^{\ast })& A(x{)}^T\\ A(x)& O\end{array}\right]\sim \left[\begin{array}{cc}{\nabla }_{xx}^2\mathcal{L}(x,{\lambda }^{\ast })& A(x{)}^T\\ O& -A(x){\left[{\nabla }_{xx}^2\mathcal{L}(x,{\lambda }^{\ast })\right]}^{-1}A(x{)}^T\end{array}\right],$$此阵的行列式不为0. 因此此系统可视作牛顿方程的良态重构(well conditioned reformulation). 不过, 由于${\nabla }^{2}f(x)+{\sum }_{i\in \mathcal{E}}{\mu }_k{c}_i(x){\nabla }^2{c}_i(x)$对${\nabla }_{xx}^2\mathcal{L}(x,{\lambda }^{\ast })$(也即${\mu }_k{c}_i(x)$对$-{\lambda }_i^{\ast }$)的近似可能不会很好, 所以这两个系统可能都不会得出一个比较好的搜索方向$p$. 这一事实可能就说明二次模型并不是$Q(\cdot ;{\mu }_k)$的一个合适的近似模型, 而牛顿步可能本身就不适合作为搜索方向. 后续我们会讨论针对以上的补救措施.

如果要通过求解重构的问题计算迭代步, 我们需要求解一个$n+|\mathcal{E}|$维而不是$n$维的系统. 后面在下一章中讨论逐步二次规划(SQP)时我们也得求解类似的系统. 事实上,

• 当${\mu }_k$很大, 上述重构可视作SQP迭代步的正则化, 其中的$-(1/{\mu }_k)I$使得迭代矩阵即使$A(x)$行亏秩时依然是非奇异的.
• 当${\mu }_k$较小, 则重构系统的第二行说明计算出的迭代步并不能较好地满足约束的线性化. 我们是不希望这种情况发生的, 因为这样一来迭代步就不能朝着可行性改进多少, 从而全局表现不够高效.
• 若 { μ k } \{\mu_k\} {μk​}趋于无穷过快, 我们就又可能丧失线性化满足时获得超线性收敛速度的机会. 具体讨论可见下一章.

总之, 重构系统可以视作无约束极小在二次惩罚函数$Q(\cdot ;{\mu }_k)$上的应用, 也可以视作SQP方法的一个变体.

2. 非光滑精确罚函数法

有些罚函数是精确的(exact), 即对于特定的惩罚因子, 仅做一次极小化即可得到非线性规划问题的精确解. 这一性质是很吸引人的, 因为它使得罚函数方法不依赖于惩罚因子的更新策略. 第1节介绍的二次罚函数并不是精确的, 这是因为对任何$\mu$的正值, 其极小点一般都不是原非线性规划的解. 本节我们讨论非光滑精确罚函数法.
对一般非线性规划的一种广受欢迎的非光滑罚函数是$l_1$罚函数, 定义为$${\varphi }_1(x;\mu )=f(x)+\mu \sum _{i\in \mathcal{E}}|c_i(x)|+\mu \sum _{i\in \mathcal{I}}[c_i(x){]}^{-}.$$其称谓源于它使用的惩罚项是$\mu$乘以约束违反度的$l_1$范数. 注意${\varphi }_1(x;\mu )$在某些$x$上并不可微.

下面的定理揭示了$l_1$罚函数的精确性.

定理3 设$x^{\ast }$为非线性规划问题的严格局部解, 且其满足一阶必要性条件, 有Lagrange乘子${\lambda }_i^{\ast },i\in \mathcal{E}\cup \mathcal{I}$. 于是$x^{\ast }$为${\varphi }_1(x;\mu ),\forall \mu >{\mu }^{\ast }$的局部极小点, 其中$${\mu }^{\ast }=\Vert {\lambda }^{\ast }{\Vert }_{\infty }=\underset{i\in \mathcal{E}\cup \mathcal{I}}{\max }|{\lambda }_i^{\ast }|.$$另外, 若在$x^{\ast }$还满足二阶充分性条件, 且$\mu >{\mu }^{\ast }$, 则$x^{\ast }$为${\varphi }_1(x;\mu )$的严格局部极小点.

证明可见 Han S P和Mangasarian O L¹的定理4.4.

粗略地说, 在非线性规划的解$x^{\ast }$处, 任何向不可行域的移动都会被严重地惩罚, 这样使得函数值会大于${\varphi }_1(x^{\ast };\mu )=f(x^{\ast })$(即方向导数都大于0), 强迫$\varphi (\cdot ;\mu )$的极小点落在$x^{\ast }$.

例2 考虑如下单变量问题:$$\min x,\mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}x\ge 1,$$其解为$x^{\ast }=1$. 我们有$${\varphi }_1(x;\mu )=x+\mu [x-1{]}^{-}=\{\begin{array}{ll}(1-\mu )x+\mu ,& x\le 1,\\ x& x>1.\end{array}$$于是从下图可见, 当$\mu >1$, 罚函数有极小点$x^{\ast }=1$, 但在$\mu <1$时却是个单调递增函数.

由于罚函数法通过直接极小化罚函数起作用, 所以我们需要刻画${\varphi }_1$的稳定点. 即使${\varphi }_1$不可微, 但它沿任何方向具有方向导数$D({\varphi }_1(x;\mu );p)$.

定义1 一点$\hat{x}\in {\mathbb{R}}^n$是罚函数${\varphi }_1(x;\mu )$的稳定点, 若$$D({\varphi }_1(\hat{x};\mu );p)\ge 0,\forall p\in {\mathbb{R}}^n.$$类似地, $\hat{x}$为不可行度量(the measure of infeasibility)$$h(x)=\sum _{i\in \mathcal{E}}|c_i(x)|+\sum _{i\in \mathcal{I}}[c_i(x){]}^{-}$$的稳定点, 若$D(h(\hat{x});p)\ge 0,\forall p\in {\mathbb{R}}^n$. 若一点不可行但对不可行度量$h$是稳定的, 则我们称其为不可行稳定点(infeasible stationary point).

对于例2中的函数, 我们在$x^{\ast }=1$有$$D({\varphi }_1(x^{\ast };\mu );p)=\{\begin{array}{ll}p,& p\ge 0,\\ (1-\mu )p,& p<0.\end{array}$$于是当$\mu >1$, $D({\varphi }_1(x^{\ast };\mu );p)\ge 0,\forall p\in \mathbb{R}$.

下面的结论则在一定程度上弥补了定理3的不足, 即证明了反方向: 在一定条件下, ${\varphi }_1(x;\mu )$的稳定点对应了原非线性规划的KKT点.

定理4 设$\hat{x}$为罚函数${\varphi }_1(x;\mu ),\forall \mu >\hat{\mu }>0$的稳定点. 若$\hat{x}$为原非线性规划的可行点, 则其为KKT点. 若不可行, 则是不可行稳定点.

证明: 若$\hat{x}$可行, 于是$$D({\varphi }_1(\hat{x};\mu );p)=\nabla f(\hat{x}{)}^{T}p+\mu \sum _{i\in \mathcal{E}}\left|\nabla c_i(\hat{x}{)}^{T}p\right|+\mu \sum _{i\in \mathcal{I}\cap \mathcal{A}(\hat{x})}{\left[\nabla c_i(\hat{x}{)}^{T}p\right]}^{-},???$$其中$\mathcal{A}(\hat{x})$为在$\hat{x}$处的积极集. 考虑任何线性化可行方向集$\mathcal{F}(\hat{x})$中的$p$, 由$\mathcal{F}(\hat{x})$的定义, 我们有$$\left|\nabla c_i(\hat{x}{)}^{T}p\right|+\sum _{i\in \mathcal{I}\cap \mathcal{A}(\hat{x})}{\left[\nabla c_i(\hat{x}{)}^{T}p\right]}^{-}=0,$$因此由稳定性假设, 我们有$$0\le D({\varphi }_1(\hat{x};\mu );p)=\nabla f(\hat{x}{)}^{T}p,\forall p\in \mathcal{F}(\hat{x}).$$利用Farkas引理, 可得存在${\hat{\lambda }}_i:{\hat{\lambda }}_i\ge 0,i\in \mathcal{I}\cap \mathcal{A}(\hat{x})$, 使得$$\nabla f(\hat{x})=\sum _{i\in \mathcal{A}(\hat{x})}{\hat{\lambda }}_i\nabla c_i(\hat{x}).$$这说明$\hat{x}$为KKT点.

若$\hat{x}$不可行, 反证, 假设存在$p$使得$D(h(\hat{x});p)<0$. 由$$0\le D({\varphi }_1(\hat{x};\mu );p)=\nabla f(\hat{x}{)}^{T}p+\mu D(h(\hat{x});p),\forall \mu 充分大,$$可推得$\nabla f(\hat{x}{)}^{T}p\to \infty$, 矛盾! 证毕.

例3 再考虑例1中的问题, 这里使用$l_1$罚函数,$${\varphi }_1(x;\mu )=x_1+x_2+\mu |x_1^2+x_2^2-2|.$$下图描绘了$\varphi (x;2)$的等高线, 其极小点可见即为原问题的解$x^{\ast }=(-1,-1{)}^T$.

事实上由定理3, 我们知道对于所有的$\mu >|{\lambda }^{\ast }|=0.5$, ${\varphi }_1(x;\mu )$的极小点都是$x^{\ast }$. 而图中等高线的尖角则体现了$l_1$罚函数的不光滑性.

这些结论为基于$l_1$罚函数的算法框架提供了动机.

框架2 (Classical $l_1$ Penalty Method)
Given ${\mu }_0>0$, tolerance $\tau >0$, starting point $x_0^s$;
for $k=0,1,2,\dots$
\quad\quad Find an approximate minimizer $x_k$ of ${\varphi }_1(x;{\mu }_k)$, starting at $x_k^s$;
\quad\quad if $h(x_k)\le \tau$
\quad\quad\quad\quad stop with approximate solution $x_k$;
\quad\quad end (if)
\quad\quad Choose new penalty parameter ${\mu }_{k+1}>{\mu }_k$;
\quad\quad Choose new starting point $x_{k+1}^s$;
end (for)

这其中对于${\varphi }_1(x;{\mu }_k)$的极小因非光滑性而变得难以实现. 但若使用${\varphi }_1(x;{\mu }_k)$的光滑模型近似, 再对光滑模型做极小, 事情就变得好理解多了. 我们将在下一小节介绍, 这有点类似于逐步二次规划算法.

更新惩罚因子${\mu }_k$的最简单格式是乘以某个常数(比如5或10). 此格式有时很好用, 但也可能会影响效率; 若初始的惩罚因子${\mu }_0$过小, 则框架2可能需要更多的迭代更新才能取到一个比较合适的因子; 在迭代早期, 迭代点可能会逐渐远离解$x^{\ast }$, 这时就需要尽早终止${\varphi }_1(x;{\mu }_k)$的极小化过程, $x_k^s$也应该重新设置为一个先前的迭代点; 若${\mu }_k$相当大, 则极小化罚函数的困难度也会增加, 需要的迭代书也会变多. 之后我们会讨论如何选取惩罚因子.

2.1 实用$l_1$罚函数法

如前所述, ${\varphi }_1(x;\mu )$在任何$x:c_i(x)=0,\exists i\in \mathcal{E}\cup \mathcal{I}$都不可微. 考虑到此函数不可微的特殊性, 我们不考虑非光滑优化的技术(比如捆集法(bundle methods)), 而是充分考虑特殊性, 重新设计算法. 在无约束优化中, 我们构建${\varphi }_1(x;\mu )$的简单模型并通过极小化模型得到迭代步. 而这一简单模型可以通过线性化约束$c_i$并以二次函数近似非线性$f$得到:$$q(p;\mu )=f(x)+\nabla f(x{)}^{T}p+\frac{1}{2}{p}^{T}Wp+\mu \sum _{i\in \mathcal{E}}|c_i(x)+\nabla c_i(x{)}^{T}p|+\mu \sum _{i\in \mathcal{I}}[c_i(x)+\nabla c_i(x{)}^{T}p{]}^{-},$$其中$W$为对称阵, 其通常包含了$f$和$c_i,i\in \mathcal{E}\cup \mathcal{I}$的二阶导数信息. 模型$q(p;\mu )$仍不光滑, 但我们可以添加人工变量, 构造光滑的二次规划问题:$$\begin{array}{rl}\underset{p,r,s,t}{\min }& f(x)+\frac{1}{2}{p}^{T}Wp+\nabla f(x{)}^{T}p+\mu \sum _{i\in \mathcal{E}}(r_i+s_i)+\mu \sum _{i\in \mathcal{I}}{t}_i\\ \mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}& \nabla c_i(x{)}^{T}p+c_i(x)=r_i-s_i,i\in \mathcal{E},\\ & \nabla c_i(x{)}^{T}p+c_i(x)\ge -t_i,i\in \mathcal{I},\\ & r,s,t\ge 0.\end{array}$$此子问题可用标准的二次规划求解器求解. 即使再添加一个盒型(box-shaped)信赖域约束$\Vert p{\Vert }_{\infty }\le \Delta$, 它仍是一个二次规划. 这样极小化${\varphi }_1$的方法与逐步二次规划联系紧密, 我们将在下一章讨论之.

选取和更新惩罚因子${\mu }_k$的策略对迭代的成功与否很关键. 我们曾提到过一种简单但不一定高效的方法, 是先选取初值, 再重复地增大它直到可行性条件满足. 在此方法的一些变体中, 每步迭代惩罚因子的选取满足${\mu }_k>\Vert {\lambda }_k{\Vert }_{\infty }$, 其中${\lambda }_k$为$x_k$处Lagrange乘子的估计. 此策略是基于定理2的, 但由于乘子估计可能不够精确, 且在离解不够近时往往无法得到较好的${\mu }_k$, 因此这种策略也并不总是能成功.

因选取合适${\mu }_k$较为困难, 非光滑罚函数法在上世纪九十年代逐渐淡出, 从而刺激了对滤子法的发展. 相比较而言, 滤子法不需要选取惩罚因子. 不过近些年来, 又出现了对罚函数法新一轮的研究兴趣, 部分由于它们能够处理退化的问题. 而新的更新惩罚因子的方法在某些特定的情形下, 能够克服先前的困难. 具体可见下一章的算法5.

在选取初始点$x_{k+1}^s$时需要尤其注意. 若${\mu }_k$对当前合适, 则我们可以设$x_{k+1}^s$为${\varphi }_1(x;{\mu }_k)$的极小点$x_k$. 否则应考虑选取更早些迭代的极小点.

2.2 一般的非光滑罚函数法

精确非光滑罚函数可以定义为除了$l_1$范数以外的其他范数形式. 我们可以写作$$\varphi (x;\mu )=f(x)+\mu \Vert c_{\mathcal{E}}(x)\Vert +\mu \Vert [c_{\mathcal{I}}(x){]}^{-}\Vert ,$$其中$\Vert \cdot \Vert$为任意向量范数, 且所有的等式、不等式约束都分别存放于向量值函数$c_{\mathcal{E}},c_{\mathcal{I}}$中. 框架2可用于任何此类罚函数. 我们仅需重新定义不可行性度量为$h(x)=\Vert c_{\mathcal{E}}(x)\Vert +\Vert [c_{\mathcal{I}}(x){]}^{-}\Vert$. 事实上可以证明, 这些范数所对应的不同罚函数的极小点是一一对应的, 证明可见Han S P和Mangasarian O L²的定理4.2. 最常用的范数为$l_1,l_{\infty },l_2$(无平方)范数. 对$l_{\infty }$范数容易得到类似于上一小节的线性化重构.

对$l_1$罚函数成立的理论性质对一般形式也成立. 在定理3中, 我们把不等式换成$$\mu \ge {\mu }^{\ast }=\Vert {\lambda }^{\ast }{\Vert }_D,$$其中$\Vert \cdot {\Vert }_D$是$\Vert \cdot \Vert$的对偶范数. 而定理4则无需修改直接可应用.

下说明此类的罚函数若是精确的必是非光滑的. 为说明简单, 我们仅考虑只有一个等式约束$c_1(x)=0$的情形, 考虑罚函数$$\varphi (x;\mu )=f(x)+\mu h(c_i(x)),$$其中$h:\mathbb{R}\to \mathbb{R}$为满足$h(y)\ge 0,\forall y\in \mathbb{R};h(0)=0$的函数. 若$h$连续可微, 则$h$有一极小点在$0$. 于是$\nabla h(0)=0$. 若$x^{\ast }$为问题的局部解, 则$c_i(x^{\ast })=0\Rightarrow \nabla h(c_1(x^{\ast }))=0$. 若$x^{\ast }$为$\varphi (x;\mu )$的局部极小点(只需$\mu$足够大), 则有$$0=\nabla \varphi (x^{\ast };\mu )=\nabla f(x^{\ast })+\mu \nabla c_i(x^{\ast })\nabla h(c_1(x^{\ast }))=\nabla f(x^{\ast }).$$然而, $f$在约束优化问题中的解一般不是稳定点, 因此$h$连续可微的假设错误, $\varphi (\cdot ;\mu )$不可能光滑.

非光滑罚函数也可用作其他机制下计算迭代步时的价值函数. 具体可见第十五章第四节的一般性讨论和第十八和第十九章中的具体实施.

3. 增广Lagrange函数法: 等式约束情形

下面我们讨论一种被称为乘子法the method of multipliers或增广Lagrange函数法augmented Lagragian methods的方法. 此法与第1节中的二次罚函数法相关, 但它通过显式引入Lagrange乘子估计至目标函数, 减小了问题的病态程度. 引入显式Lagrange乘子估计后的函数就被称为增广Lagrange函数(augmented Lagrangian function). 与第2节中讨论的罚函数相反, 增广Lagrange函数极大地保留了光滑性, 因此具体实施可调用标准的无约束或带界约束的优化软件求解.
本节我们在Lagrange乘子估计上, 用上标表示迭代指标, 下标表示分量. 而对于其他的变量, 使用下标表示迭代指标.

3.1 动机和算法框架

我们首先考虑等式约束问题. 二次罚函数$Q(x;\mu )$以不可行度量的平方惩罚违反行为. 但我们从定理2中也看到了, $Q(x;{\mu }_k)$的极小点$x_k$并不满足可行性条件$c_i(x)=0,i\in \mathcal{E}$, 而是近似地有$$c_i(x_k)\approx -{\lambda }_i^{\ast }/{\mu }_k,\forall i\in \mathcal{E}.$$当然, $c_i(x_k)\to 0,{\mu }_k\to \infty$, 但考虑到条件数的问题, 我们可能会考虑能否修改函数$Q(x;{\mu }_k)$的形式来避免这种系统上的偏差, 即使得近似极小点无需${\mu }_k$充分大就差不多满足等式约束$c_i(x)=0$.

增广Lagrange函数${\mathcal{L}}_A(x,\lambda ;\mu )$正满足我们的需求. 基于定理2, 它向目标加入了Lagrange乘子$\lambda$的显式估计, 即有$${\mathcal{L}}_A(x,\lambda ;\mu )\stackrel{def}{=}f(x)-\sum _{i\in \mathcal{E}}{\lambda }_i{c}_i(x)+\frac{\mu }{2}\sum _{i\in \mathcal{E}}{c}_i^2(x).$$以上形式, 与标准的Lagrange函数相比多了平方项, 而与二次罚函数相比则多了关于$\lambda$的求和项. 从这个角度上看, 它是Lagrange函数和二次罚函数的结合.

下面我们设计一种对$x$极小化的算法, 其中固定惩罚因子$\mu$为${\mu }_k>0$, 固定$\lambda$为${\lambda }^k$. 用$x_k$表示${\mathcal{L}}_A(x,{\lambda }^k;{\mu }_k)$的近似极小点, 于是由无约束极小的最优性条件可知$$0={\nabla }_x{\mathcal{L}}_A(x_k,{\lambda }^k;{\mu }_k)=\nabla f(x_k)-\sum _{i\in \mathcal{E}}[{\lambda }_i^k-{\mu }_k{c}_i(x_k)]\nabla c_i(x_k).$$比较约束优化的一阶最优性条件, 若$\nabla c_i(x_k)$线性无关, 则有$${\lambda }_i^{\ast }\approx {\lambda }_i^k-{\mu }_k{c}_i(x_k),\forall i\in \mathcal{E}.$$整理可得$$c_i(x_k)\approx -\frac{1}{{\mu }_k}({\lambda }_i^{\ast }-{\lambda }_i^k),\forall i\in \mathcal{E},$$从这可看出, 若${\lambda }^k$接近于最优乘子${\lambda }^{\ast }$, 则$x_k$处的不可行性度量将会小于$(1/{\mu }_k)$, 而不是成比例(见定理2). 而上面的式子也为我们提供了更新当前估计${\lambda }^k$的公式:$${\lambda }_i^{k+1}={\lambda }_i^k-{\mu }_k{c}_i(x_k),\forall i\in \mathcal{E}.$$

于是我们得到下述算法框架.

框架3 (Augmented Lagrangian Method-Equality Constraints)
Given ${\mu }_0>0$, tolerance ${\tau }_0>0$, starting points $x_0^s$ and ${\lambda }^0$;
for $k=0,1,2,\dots$
\quad\quad Find an approximate minimizer $x_k$ of ${\mathcal{L}}_A(\cdot ,{\lambda }^k;{\mu }_k)$, starting at $x_k^s$, and terminating when $\Vert {\nabla }_x{\mathcal{L}}_A(x_k,{\lambda }^k;{\mu }_k)\Vert \le {\tau }_k$;
\quad\quad if a convergence test for the original problem is satisfied
\quad\quad\quad\quad stop with approximate solution $x_k$;
\quad\quad end (if)
\quad\quad Update Lagrange multipliers to obtain ${\lambda }^{k+1}$;
\quad\quad Choose new penalty perameter ${\mu }_{k+1}\ge {\mu }_k$;
\quad\quad Set starting point for the next iteration to $x_{k+1}^s=x_k$;
\quad\quad Select tolerance ${\tau }_{k+1}$;
end (for)

我们通过一个例子来阐释此法可以无需将$\mu$增加到$\infty$即可保证收敛. 于是病态的情况要比框架1有所缓解; 初始点$x_{k+1}^s$的选取也就不那么关键; 容忍限${\tau }_k$可依据不可行性${\sum }_{i\in \mathcal{E}}|c(x_k)|$选取; 若当前迭代下不可行性度量没有充分下降, 则惩罚因子$\mu$可能就需要增大一些.

例4 再次考虑例1中的问题, 其增广Lagrange函数为$${\mathcal{L}}_A(x,\lambda ;\mu )=x_1+x_2-\lambda (x_1^2+x_2^2-2)+\frac{\mu }{2}(x_1^2+x_2^2-2{)}^2.$$如前所述, 原问题解为$x^{\ast }=(-1,-1{)}^T$, 最优Lagrange乘子为${\lambda }^{\ast }=-0.5$.
假设在第$k$步迭代, ${\mu }_k=1$, ${\lambda }^k=-0.4$. 增广Lagrange函数的等高线可见下图. 与之对比, ${\mu }_k=1$的二次罚函数等高线图可见图1.

注意到等高线之间的跨度说明此问题的条件数与二次罚函数$Q(x;1)$差不多. 但此时的极小点$x_k\approx (-1.02,-1.02{)}^T$则更靠近解$x^{\ast }=(-1,-1{)}^T$. 这个例子表明, 将Lagrange乘子估计项包含进${\mathcal{L}}_A(x,\lambda ;\mu )$中可以得到相较于二次罚函数法的精度提升.

3.2 增广Lagrange函数的性质

下面我们证明关于增广Lagrange函数使用和等式约束问题乘子法的两条结论.

第一条证实, 若我们有精确Lagrange乘子${\lambda }^{\ast }$的信息, 则原问题的解$x^{\ast }$是${\mathcal{L}}_A(x,{\lambda }^{\ast };\mu ),\mu$充分大的严格极小点. 尽管我们一般不能预先知道准确的${\lambda }^{\ast }$, 但此结论和其证明表明只要$\lambda$是${\lambda }^{\ast }$的一个较好的估计, 我们可以通过在$\mu$不是特别大时, 极小化${\mathcal{L}}_A(x,\lambda ;\mu )$得到关于$x^{\ast }$的一个好估计.

定理5 令$x^{\ast }$为原问题的局部解, 且LICQ成立, 对$\lambda ={\lambda }^{\ast }$有二阶充分性条件成立. 则存在$\bar{\mu }$使得对$\forall \mu \ge \bar{\mu }$, $x^{\ast }$为${\mathcal{L}}_A(x,{\lambda }^{\ast };\mu )$的严格局部极小点.

证明: 我们证明$x^{\ast }$对$\mu$充分大, 满足关于${\mathcal{L}}_A(x,{\lambda }^{\ast };\mu )$的二阶充分性条件, 即$${\nabla }_x{\mathcal{L}}_A(x^{\ast },{\lambda }^{\ast };\mu )=0,{\nabla }_{xx}^2{\mathcal{L}}_A(x^{\ast },{\lambda }^{\ast };\mu )正定.$$由于$x^{\ast }$为原问题局部解, 且LICQ成立, 于是由KKT条件知${\nabla }_x\mathcal{L}(x^{\ast },{\lambda }^{\ast })=0,c_i(x^{\ast })=0,\forall i\in \mathcal{E}$, 从而$$\begin{array}{rl}{\nabla }_x{\mathcal{L}}_A(x^{\ast },{\lambda }^{\ast };\mu )& =\nabla f(x^{\ast })-\sum _{i\in \mathcal{E}}[{\lambda }_i^{\ast }-\mu c_i(x^{\ast })]\nabla c_i(x^{\ast })\\ & =\nabla f(x^{\ast })-\sum _{i\in \mathcal{E}}{\lambda }_i^{\ast }\nabla c_i(x^{\ast })={\nabla }_x\mathcal{L}(x^{\ast },{\lambda }^{\ast })=0,\end{array}$$这一点与$\mu$无关. 对于二阶充分性条件, 定义$A$为约束梯度矩阵, 于是$${\nabla }_{xx}^2{\mathcal{L}}_A(x^{\ast },{\lambda }^{\ast };\mu )={\nabla }_{xx}^2\mathcal{L}(x^{\ast },{\lambda }^{\ast })+\mu A^{T}A.$$若对于充分大$\mu$二阶充分性条件不成立, 则对每个$k\ge 1$, 我们可选取向量$w_k:\Vert w_k\Vert =1$使得$$0\ge w_k^T{\nabla }_{xx}^2{\mathcal{L}}_A(x^{\ast },{\lambda }^{\ast };k)w_k=w_k^T{\nabla }_{xx}^2\mathcal{L}(x^{\ast },{\lambda }^{\ast })w_k+k\Vert A{w}_k{\Vert }_2^2,$$因此$$\Vert A{w}_k{\Vert }_2^2\le -(1/k)w_k^T{\nabla }_{xx}^2\mathcal{L}(x^{\ast },{\lambda }^{\ast })w_k\to 0,k\to \infty .$$由于 { w k } \{w_k\} {wk​}位于一紧集中, 于是存在聚点$w$. 而上面的式子说明$Aw=0$. 这说明$w\in \mathcal{F}(x^{\ast })$. 进一步, $$w_k^T{\nabla }_{xx}^2\mathcal{L}(x^{\ast },{\lambda }^{\ast })w_k\le -k\Vert A{w}_k{\Vert }_2^2\le 0,$$取极限有$w^T{\nabla }_{xx}^2\mathcal{L}(x^{\ast },{\lambda }^{\ast })w\le 0$. 这与原问题的二阶充分性矛盾. 证毕.

第二条结论则立意于更实际的情形$\lambda \ne {\lambda }^{\ast }$. 它给出了${\mathcal{L}}_A(x,\lambda ;\mu )$极小点靠近$x^{\ast }$的条件, 且给出了$x_k$和更新的乘子估计${\lambda }^{k+1}$分别与$x^{\ast },{\lambda }^{\ast }$的误差界.

定理6 设定理5的条件在$x^{\ast },{\lambda }^{\ast }$上成立, 且令$\bar{\mu }$为定理5中的阈值. 则存在正标量$\delta ,ϵ,M$使得以下断言成立:

1. 对所有满足$$\Vert {\lambda }^k-{\lambda }^{\ast }\Vert \le {\mu }_k\delta ,{\mu }_k\ge \bar{\mu }$$的${\lambda }^k,{\mu }_k$, 问题$$\underset{x}{\min }{\mathcal{L}}_A(x,{\lambda }^k;{\mu }_k),\mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}\Vert x-x^{\ast }\Vert \le ϵ$$有唯一解$x_k$. 进一步地, $$\Vert x_k-x^{\ast }\Vert \le M\Vert {\lambda }^k-{\lambda }^{\ast }\Vert /{\mu }_k.$$
2. 对所有满足$$\Vert {\lambda }^k-{\lambda }^{\ast }\Vert \le {\mu }_k\delta ,{\mu }_k\ge \bar{\mu }$$的${\lambda }^k,{\mu }_k$, 我们有$$\Vert {\lambda }^{k+1}-{\lambda }^{\ast }\Vert \le M\Vert {\lambda }^k-{\lambda }^{\ast }\Vert /{\mu }_k,$$其中${\lambda }^{k+1}$由公式$${\lambda }_i^{k+1}={\lambda }_i^k-{\mu }_k{c}_i(x_k),\forall i\in \mathcal{E}$$给出.
3. 对所有满足$$\Vert {\lambda }^k-{\lambda }^{\ast }\Vert \le {\mu }_k\delta ,{\mu }_k\ge \bar{\mu }$$的${\lambda }^k,{\mu }_k$, 矩阵${\nabla }_{xx}^2{\mathcal{L}}_A(x_k,{\lambda }^k;{\mu }_k)$正定, 约束梯度$\nabla c_i(x_k),i\in \mathcal{E}$线性无关.

证明: 以下以$\mu ,\lambda ,\tilde{\lambda }$分别表示${\mu }_k,{\lambda }^k,{\lambda }^{k+1}$. 对$\mu >0$, 考虑以下以$(x,\tilde{\lambda },\lambda ,\mu )$为变量的系统:$$\nabla f(x)-\nabla h(x)\tilde{\lambda }=0,h(x)+(\lambda -\tilde{\lambda })/\mu =0,$$其中$$h(x)=\left[\begin{array}{c}{c}_1(x)\\ c_2(x)\\ ⋮\\ c_m(x)\end{array}\right],\nabla h(x)=\left[\begin{array}{cccc}\nabla c_1(x)& \nabla c_2(x)& \cdots & \nabla c_m(x)\end{array}\right].$$引入变量$t\in {\mathbb{R}}^m,\gamma \in \mathbb{R}$:$$t=(\lambda -{\lambda }^{\ast })/\mu ,\gamma =1/\mu ,$$上述系统可写作$$\nabla f(x)-\nabla h(x)\tilde{\lambda }=0,h(x)+t+\gamma {\lambda }^{\ast }-\gamma \tilde{\lambda }=0.$$对$t=0,\gamma \in [0,1/\bar{\mu }]$, 上述系统有解$x=x^{\ast },\tilde{\lambda }={\lambda }^{\ast }$. 在此对变量$(x,\tilde{\lambda })$的Jacobi矩阵为$$\left[\begin{array}{cc}{\nabla }_{xx}^2\mathcal{L}(x^{\ast },{\lambda }^{\ast })& \nabla h(x^{\ast })\\ -\nabla h(x^{\ast }{)}^T& -\gamma I\end{array}\right].$$我们证明, 这一矩阵对所有$\gamma \in [0,1/\bar{\mu }]$可逆. 显然$\gamma =0$时这是成立的. 为证明这对于$\forall \gamma \in (0,1/\bar{\mu }]$亦成立, 设对某个$z\in {\mathbb{R}}^n,w\in {\mathbb{R}}^m$, 我们有$$\left[\begin{array}{cc}{\nabla }_{xx}^2\mathcal{L}(x^{\ast },{\lambda }^{\ast })& \nabla h(x^{\ast })\\ -\nabla h(x^{\ast }{)}^T& -\gamma I\end{array}\right]\left[\begin{array}{c}z\\ w\end{array}\right]=0$$或等价地, $$\begin{array}{rl}{\nabla }_{xx}^2\mathcal{L}(x^{\ast },{\lambda }^{\ast })z+\nabla h(x^{\ast })w& =0,\\ -\nabla h(x^{\ast }{)}^{T}z-\gamma w& =0.\end{array}$$消去$w$, 得到$$[{\nabla }_{xx}^2\mathcal{L}(x^{\ast },{\lambda }^{\ast })-(1/\gamma )\nabla h(x^{\ast })\nabla h(x^{\ast }{)}^T]z=0.$$对于$\gamma =1/\mu ,\mu \ge \bar{\mu }$, 这就得出${\nabla }_{xx}^2{\mathcal{L}}_A(x^{\ast },{\lambda }^{\ast };\mu )z=0$, 而由定理5知${\nabla }_{xx}^2{\mathcal{L}}_A(x^{\ast },{\lambda }^{\ast };\mu )>0,\mu \ge \bar{\mu }$, 所以$z=0\Rightarrow w=0$. 这就说明矩阵对所有$\gamma \in [0,1/\bar{\mu }]$可逆.
由隐函数定理, 存在$ϵ>0,\delta >0$和定义在$S(K;\delta )$(其中 K = { ( 0 , γ ) ∣ γ ∈ [ 0 , 1 / μ ˉ ] } K=\{(0,\gamma)\mid\gamma\in[0,1/\bar\mu]\} K={(0,γ)∣γ∈[0,1/μˉ​]})上的连续可微函数$\hat{x}(t,\gamma ),\hat{\lambda }(t,\gamma )$使得$${\left(|\hat{x}(t,\gamma )-x^{\ast }{|}^2+|\hat{\lambda }(t,\gamma )-{\lambda }^{\ast }{|}^2\right)}^{1/2}<ϵ,\forall (t,\gamma )\in S(K;\delta ),$$且满足$$\begin{array}{rl}\nabla f[\hat{x}(t,\gamma )]-\nabla h[\hat{x}(t,\gamma )]\hat{\lambda }(t,\gamma )& =0,\\ h[\hat{x}(t,\gamma )]+t+\gamma {\lambda }^{\ast }-\gamma \hat{\lambda }(t,\gamma )& =0.\end{array}$$显然$\delta ,ϵ$的选取可使$\nabla h[\hat{x}(t,\gamma )]$秩为$m$且$${\nabla }_{xx}^2\mathcal{L}[\hat{x}(t,\gamma ),\hat{\lambda }(t,\gamma )]-\mu \nabla h[\hat{x}(t,\gamma )]\nabla h[\hat{x}(t,\gamma ){]}^T>0,\forall (t,\gamma )\in S(K;\delta ),\mu \ge \bar{\mu }.$$对$\mu \ge \bar{\mu },\Vert \lambda -{\lambda }^{\ast }\Vert \le \delta \mu$, 定义$$x(\lambda ,\mu )=\hat{x}\left(\frac{\lambda -{\lambda }^{\ast }}{\mu },\frac{1}{\mu }\right),\tilde{\lambda }(\lambda ,\mu )=\hat{\lambda }\left(\frac{\lambda -{\lambda }^{\ast }}{\mu },\frac{1}{\mu }\right).$$于是对所有满足$$\Vert \lambda -{\lambda }^{\ast }\Vert \le \mu \delta ,\mu \ge \bar{\mu }$$的${\lambda }^k,{\mu }_k$, 有$$\begin{array}{rl}\nabla f[x(\lambda ,\mu )]-\nabla h[x(\lambda ,c)]\tilde{\lambda }(\lambda ,\mu )& =0,\\ \tilde{\lambda }(\lambda ,\mu )& =\lambda -\mu h[x(\lambda ,\mu )],\\ {\nabla }_{xx}^2\mathcal{L}[x(\lambda ,\mu ),\tilde{\lambda }(\lambda ,\mu )]-\mu \nabla h[x(\lambda ,\mu )]\nabla h[x(\lambda ,\mu ){]}^T={\nabla }_{xx}^2{\mathcal{L}}_A[x(\lambda ,\mu ),\lambda ]& >0.\end{array}$$这就证明了第三条. 而为了证明第一条和第二条, 对$$\begin{array}{rl}\nabla f[\hat{x}(t,\gamma )]-\nabla h[\hat{x}(t,\gamma )]\hat{\lambda }(t,\gamma )& =0,\\ h[\hat{x}(t,\gamma )]+t+\gamma {\lambda }^{\ast }-\gamma \hat{\lambda }(t,\gamma )& =0.\end{array}$$中的$t,\gamma$求导, 得到$$\left[\begin{array}{cc}{\nabla }_t\hat{x}(t,\gamma {)}^T& {\nabla }_{\gamma }\hat{x}(t,\gamma {)}^T\\ {\nabla }_t\hat{\lambda }(t,\gamma {)}^T& {\nabla }_{\gamma }\hat{\lambda }(t,\gamma {)}^T\end{array}\right]=A(t,\gamma )\left[\begin{array}{cc}O& O\\ -I& \hat{\lambda }(t,\gamma )-{\lambda }^{\ast }\end{array}\right],$$其中$$A(t,\gamma )={\left[\begin{array}{cc}{\nabla }_{xx}^2\mathcal{L}[\hat{x}(t,\gamma ),\hat{\lambda }(t,\gamma )]& -\nabla h[\hat{x}(t,\gamma )]\\ \nabla h[\hat{x}(t,\gamma ){]}^T& -\gamma I\end{array}\right]}^{-1}.$$对$\forall (t,\gamma ):|t|<\delta ,\gamma \in [0,1/\bar{\mu }]$, 我们有$$\begin{array}{rl}\left[\begin{array}{c}\hat{x}(t,\gamma )-x^{\ast }\\ \hat{\lambda }(t,\gamma )-{\lambda }^{\ast }\end{array}\right]& =\left[\begin{array}{c}\hat{x}(t,\gamma )-\hat{x}(0,0)\\ \hat{\lambda }(t,\gamma )-\hat{\lambda }(0,0)\end{array}\right]\\ & ={\int }_0^{1}A(\zeta t,\zeta \gamma )\left[\begin{array}{cc}O& O\\ -I& \hat{\lambda }(\zeta t,\zeta \gamma )-{\lambda }^{\ast }\end{array}\right]\left[\begin{array}{c}t\\ \gamma \end{array}\right]\mathrm{d}\zeta .\end{array}$$由前面证得的矩阵对所有$\gamma \in [0,1/\bar{\mu }]$可逆, 于是对充分小的$\delta$, $A(t,\gamma )$在 { ( t , γ ) ∣ ∣ t ∣ < δ , γ ∈ [ 0 , 1 / μ ˉ ] } \{(t,\gamma)\mid|t|<\delta,\gamma\in[0,1/\bar\mu]\} {(t,γ)∣∣t∣<δ,γ∈[0,1/μˉ​]}上已知有界. 令$\bar{M}$为界, 即$\Vert A(t,\gamma )\Vert \le \bar{M},\forall |t|<\delta ,\gamma \in [0,1/\bar{\mu }]$, 且进一步可取$\delta$充分小以保证$\bar{M}\delta <1$. 于是从上式得到$$\begin{array}{rlr}& {\left(|\hat{x}(t,\gamma )-x^{\ast }{|}^2+|\hat{\lambda }(t,\gamma )-{\lambda }^{\ast }{|}^2\right)}^{1/2}& \le \bar{M}\left(|t|+\underset{0\le \zeta \le 1}{\max }|\hat{\lambda }(\zeta t,\zeta \gamma )-{\lambda }^{\ast }|\gamma \right).\end{array}$$由此, 对$\forall (t,\gamma ):|t|<\delta ,\gamma \in [0,1/\bar{\mu }],\gamma <\delta$,$$|\hat{\lambda }(t,\gamma )-{\lambda }^{\ast }|\le \bar{M}|t|+\bar{M}\gamma \underset{0\le \zeta \le 1}{\max }|\hat{\lambda }(\zeta t,\zeta \gamma )-{\lambda }^{\ast }|.$$将上式中的$t,\gamma$换成$\zeta t,\zeta \gamma ,\zeta \in [0,1]$, 得到$$\underset{0\le \zeta \le 1}{\max }|\hat{\lambda }(\zeta t,\zeta \gamma )-{\lambda }^{\ast })|\le \frac{\bar{M}}{1-\bar{M}\gamma }|t|.$$组合起来就有, 对$\forall (t,\gamma ):|t|<\delta ,\gamma \in [0,1/\bar{\mu }],\gamma <\delta$,$${\left(|\hat{x}(t,\gamma )-x^{\ast }{|}^2+|\hat{\lambda }(t,,\gamma )-{\lambda }^{\ast }{|}^2\right)}^{1/2}\le \left(\mu +\frac{{\mu }^2\gamma }{1-\mu \gamma }\right)|t|\le \frac{\bar{M}}{1-\bar{M}\delta }|t|.$$令$\delta$充分小, 我们有$${\left(|\hat{x}(t,\gamma )-x^{\ast }{|}^2+|\hat{\lambda }(t,,\gamma )-{\lambda }^{\ast }{|}^2\right)}^{1/2}\le 2\bar{M}|t|.$$由$t$的定义, 并令$x(\lambda ,\mu )=\hat{x}(t,\gamma ),\tilde{\lambda }(\lambda ,\mu )=\hat{\lambda }(t,\gamma )$, 我们有对 ∀ ( λ , μ ) : ∣ λ − λ ∗ ∣ < δ μ , μ > max ⁡ { μ ˉ , 1 / δ } \forall(\lambda,\mu):|\lambda-\lambda^*|<\delta\mu,\mu>\max\{\bar\mu,1/\delta\} ∀(λ,μ):∣λ−λ∗∣<δμ,μ>max{μˉ​,1/δ}, 有$$|x(\lambda ,\mu )-x^{\ast }|\le 2\bar{M}|\lambda -{\lambda }^{\ast }|/\mu ,|\tilde{\lambda }(\lambda ,\mu )-{\lambda }^{\ast }|\le 2\bar{M}|\lambda -{\lambda }^{\ast }|/\mu .$$令$M=2\bar{M}$即可得第一条和第二条对 ∀ ( λ , μ ) : ∣ λ − λ ∗ ∣ < δ μ , μ > max ⁡ { μ ˉ , 1 / δ } \forall(\lambda,\mu):|\lambda-\lambda^*|<\delta\mu,\mu>\max\{\bar\mu,1/\delta\} ∀(λ,μ):∣λ−λ∗∣<δμ,μ>max{μˉ​,1/δ}是成立的. 而由于$x(\cdot ,\cdot )$连续可微, 因此我们也可找到$M$使得第一条和第二条对 ∀ ( λ , μ ) : ∣ λ − λ ∗ ∣ < δ μ , ≤ ˉ μ ≤ max ⁡ { μ ˉ , 1 / δ } \forall(\lambda,\mu):|\lambda-\lambda^*|<\delta\mu,\bar\le\mu\le\max\{\bar\mu,1/\delta\} ∀(λ,μ):∣λ−λ∗∣<δμ,≤ˉ​μ≤max{μˉ​,1/δ}成立. 证毕.

此定理解释了增广Lagrange函数法的一些显著的性质.

1. $$\Vert x_k-x^{\ast }\Vert \le M\Vert {\lambda }^k-{\lambda }^{\ast }\Vert /{\mu }_k$$表示若${\lambda }^k$精确或惩罚因子${\mu }_k$很大, 则$x_k$就离得$x^{\ast }$很近. 因此这就给出了两种改进$x_k$精确度的方法, 而二次罚函数法只给出了一种方法——增大惩罚因子${\mu }_k$.
2. $$\Vert {\lambda }^{k+1}-{\lambda }^{\ast }\Vert \le M\Vert {\lambda }^k-{\lambda }^{\ast }\Vert /{\mu }_k$$则表明, 局部上我们可以通过选取充分大的${\mu }_k$保证乘子上精确度的提升.
3. 无约束极小的二阶充分条件在第$k$个子问题下依然是成立的, 因此使用标准的无约束极小技术是有望得到好的结果的.

4. 实用增广Lagrange函数法

本节我们讨论实用的增广Lagrange函数法, 特别地, 是要处理带不等式约束的情形. 我们分别讨论三种方式——界约束重构、约束线性化重构和无约束重构.

4.1 界约束重构

给定一般的非线性规划问题, 我们可以增加松弛变量$s_i$将问题转化为带等式约束的问题, 即$$c_i(x)-s_i=0,s_i\ge 0,\forall i\in \mathcal{I}.$$而界约束$l\le x\le u$则无需变换. 这样, 我们就可以将非线性规划写作$$\underset{x\in {\mathbb{R}}^n}{\min }f(x),\mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}{c}_i(x)=0,i=1,2,\dots ,m,l\le x\le u.$$其中下界向量$l$的一些分量可能为$-\infty$; $u$同理.

界约束重构的增广Lagrange函数法the bound-constrained Lagrangian methods, BCL methods只包含等式约束, 即$${\mathcal{L}}_A(x,\lambda ;\mu )=f(x)-\sum _{i=1}^m{\lambda }_i{c}_i(x)+\frac{\mu }{2}\sum _{i=1}^m{c}_i^2(x).$$而界约束则显式地加在子问题上, 于是有子问题$$\underset{x}{\min }{\mathcal{L}}_A(x,\lambda ;\mu ),\mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}\le x\le u.$$在此问题被近似求解后, 我们就更新乘子$\lambda$和惩罚因子$\mu$, 之后重复.

一种求解带界约束重构非线性规划的技术是非线性投影梯度法, 可见下一章. 考虑子问题的KKT条件, 我们会发现$x$为子问题解的一阶必要条件是$$x-P(x-{\nabla }_x{\mathcal{L}}_A(x,\lambda ;\mu ),l,u)=0,$$其中$P(g,l,u)$为将$g\in {\mathbb{R}}^n$投影至长方体$[l,u]$上的投影算子, 其定义为$$P(g,l,u{)}_i=\{\begin{array}{ll}{l}_i,& g_i\le l_i,\\ g_i,& g_i\in (l_i,u_i),\\ u_i,& g_i\ge u_i,\end{array}i=1,2,\dots ,n.$$以下为算法表述.

算法4 (Bound-Constrained Lagrangian Method)
Choose an initial point $x_0$ and initial multipliers ${\lambda }^0$;
Choose convergence tolerances ${\eta }_{\ast }$ and $w_{\ast }$;
Set ${\mu }_0=10,w_0=1{\mu }_0,{\eta }_0=1/{\mu }_0^{0.1}$;
for $k=0,1,2,\dots$
\quad\quad Find an approximate solution $x_k$ of the subproblem such that$$\Vert x_k-P(x_k-{\nabla }_x{\mathcal{L}}_A(x_k,{\lambda }^k;{\mu }_k),l,u)\Vert \le w_k;$$ \quad\quad if $\Vert c(x_k)\Vert \le {\eta }_k$
\quad\quad\quad\quad (test for convergence)
\quad\quad\quad\quad if $\Vert c(x_k)\Vert \le {\eta }_{\ast }$ and $\Vert x_k-P(x_k-{\nabla }_x{\mathcal{L}}_A(x_k,{\lambda }^k;{\mu }_k),l,u)\Vert \le w_{\ast }$
\quad\quad\quad\quad\quad\quad stop with approximate solution $x_k$;
\quad\quad\quad\quad end (if)
\quad\quad\quad\quad (update multipliers, tighten tolerances)
\quad\quad\quad\quad ${\lambda }^{k+1}={\lambda }^k-{\mu }_{k}c(x_k)$;
\quad\quad\quad\quad ${\mu }_{k+1}={\mu }_k$;
\quad\quad\quad\quad ${\eta }_{k+1}={\eta }_k/{\mu }_{k+1}^{0.9}$;
\quad\quad\quad\quad $w_{k+1}=w_k/{\mu }_{k+1}$;
\quad\quad else
\quad\quad\quad\quad (increase penalty parameter, tighten tolerances)
\quad\quad\quad\quad ${\lambda }^{k+1}={\lambda }^k$;
\quad\quad\quad\quad ${\mu }_{k+1}=100{\mu }_k$;
\quad\quad\quad\quad ${\eta }_{k+1}=1/{\mu }_{k+1}^{0.1}$;
\quad\quad\quad\quad $w_{k+1}=1/{\mu }_{k+1}$;
\quad\quad end (if)
end (for)

算法主要分岔位于近似求解子问题后, 算法将验证约束违反是否充分下降, 即验证条件$$\Vert c(x_k)\Vert \le {\eta }_k.$$若此条件成立, 则不改变惩罚因子, 而Lagrange乘子估计则按公式更新, 在下一次迭代前减小$w_k,{\eta }_k$. 若条件不成立, 则增大惩罚因子保证下一次迭代求解子问题时会更加“关照”约束违反度. 此时就不更新Lagrange乘子, 毕竟主要任务是改进可行性.

算法4中出现的常数0.1,0.9,100在某种程度上可以任取. 软件LANCELOT使用带信赖域的投影梯度法求解界约束非线性子问题. 此时, 投影梯度法构建增广Lagrange函数${\mathcal{L}}_A$的一个二次模型, 并通过近似求解以下信赖域问题得到迭代步$d$:$$\begin{array}{rl}\underset{d}{\min }& \frac{1}{2}{d}^T\left[{\nabla }_{xx}^2\mathcal{L}(x_k,{\lambda }^k)+{\mu }_k{A}_k^T{A}_k\right]d+{\nabla }_x{\mathcal{L}}_A(x_k,{\lambda }^k;{\mu }_k{)}^{T}d\\ \mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}& l\le x_k+d\le u,\Vert d{\Vert }_{\infty }\le \Delta ,\end{array}$$其中$A_k=A(x_k)$, $\Delta$为信赖域半径. 我们可以将信赖域约束写作界约束$-\Delta e\le d\le \Delta e$. 求解此子问题的算法每步迭代由两个阶段组成:

• 第一阶段, 做投影梯度线搜索以决定$d$的哪些分量需设为它们的边界值.
• 第二阶段, 使用CG极小化子问题, 其中不再牵涉第一阶段中被固定的分量.

这有点类似于第十六章二次规划的投影梯度法. 第二阶段相当于在可行多面体的某个面上做的极小化. 这一算法的好处在于, 无需计算KKT矩阵或者约束Jacobi矩阵$A_k$的分解. CG迭代仅需计算矩阵-向量乘积.

需要说明的是, Lagrange函数的Hessian阵${\nabla }_{xx}^2\mathcal{L}(x_k,{\lambda }^k)$可以用基于BFGS或SR1更新公式的拟牛顿近似矩阵替代. LANCELOT包也被设计成可充分利用目标函数与约束中的部分可分结构, 以高效计算Lagrange函数的Hessian阵或拟牛顿近似矩阵, 可见第七章.

4.2 约束线性化重构

约束线性化重构的增广Lagrange函数法linearly constrained Lagragian methods, LCL methods的主要想法是, 在线性化原约束的限制下, 极小化Lagrange函数(或增广Lagrange函数), 产生迭代步. 对于上一小节中的一般非线性规划$$\underset{x\in {\mathbb{R}}^n}{\min }f(x),\mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}{c}_i(x)=0,i=1,2,\dots ,m,l\le x\le u,$$LCL方法的子问题就具有以下形式:$$\begin{array}{rl}\underset{x}{\min }& F_k(x)\\ \mathrm{s}\mathrm{u}\mathrm{b}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{t}\mathrm{o}& c(x_k)+A_k(x-x_k)=0,l\le x\le u.\end{array}$$对$F_k(x)$我们有许多不同的选择. 早期的LCL方法定义$$F_k(x)=f(x)-\sum _{i=1}^m{\lambda }_i^k{\bar{c}}_i^k(x),$$其中${\lambda }^k$为当前的Lagrange乘子估计, ${\bar{c}}_i^k(x)$为$c_i(x)$和它在$x_k$处线性化的差函数, 即$${\bar{c}}_i^k(x)=c_i(x)-c_i(x_k)-\nabla c_i(x_k{)}^T(x-x_k).$$可以证明, 随着$x_k$收敛于解$x^{\ast }$, 对应于子问题中等式约束的Lagrange乘子将收敛于最优乘子. 因此, 我们可设${\lambda }^k$为前一次迭代中对应于等式约束的Lagrange乘子.

现在LCL方法定义$F_k$为增广Lagrange函数$$F_k(x)=f(x)-\sum _{i=1}^m{\lambda }_i^k{\bar{c}}_i^k(x)+\frac{\mu }{2}\sum _{i=1}^m[{\bar{c}}_i^k(x){]}^2.$$这一定义与之前相比, 在实际中具有更好的全局收敛效果.

可以看到, 上述$F_k$与之前定义的增广Lagrange函数具有很大的相似度, 而不同在于原本的约束$c_i(x)$被替换为了${\bar{c}}_i^k(x)$, 后者仅捕捉$c_i$二阶及以上的信息. 线性化约束重构子问题与原本的增广Lagrange函数子问题的不同则在于, 前者要求新迭代点$x$精确满足等式约束的线性化, 同时又将约束的线性部分从目标函数中抽出(这就是用${\bar{c}}_i^k$换$c_i$). 类似于算法4的程序可用于更新惩罚因子$\mu$以及调整容忍限.

因${\bar{c}}_i^k(x)$在$x=x_k$处梯度为0, 于是$\nabla F_k(x_k)=\nabla f(x_k)$. 我们也可验证$F_k$的Hessian与Lagrange函数或增广Lagrange函数的Hessian也是紧密相关的. 基于这些性质, 线性化重构子问题就类似于第十八章中的SQP子问题.

MINOS包使用增广形式的模型函数, 用既约投影梯度法求解线性化子问题, 其中对$F_k$的既约Hessian用拟牛顿近似. 为确保等式约束Lagrange乘子的精度, MINOS中每一次子问题的求解需要更多的目标函数值和约束函数(以及它们的梯度)值. 这要多于SQP算法或内点法. 不过, 这样一来所需求解的子问题数一般比其他方法要少一些.

4.3 无约束化重构

利用近邻点proximal point方法, 我们可以得到带不等式约束问题的无约束化重构增广Lagrange函数子问题. 为说明简便, 假设问题无等式约束, 于是我们可以将原问题等价地写作无约束优化问题$$\underset{x\in {\mathbb{R}}^n}{\min }F(x),$$其中$$F(x)=\underset{\lambda \ge 0}{\max }\left\{f(x)-\sum _{i\in \mathcal{I}}{\lambda }_i{c}_i(x)\right\}=\{\begin{array}{ll}f(x),& xisfeasible,\\ \infty ,& otherwise.\end{array}$$由于$F$不光滑, 直接极小化$F$并不实际. 于是再次想到用光滑的近似模型. 我们以$\hat{F}(x;{\lambda }^k,{\mu }_k)$代替$F$, 前者依赖于惩罚因子${\mu }_k$和Lagrange乘子估计${\lambda }^k$. 其定义为$$\hat{F}(x;{\lambda }^k,{\mu }_k)=\underset{\lambda \ge 0}{\max }\left\{f(x)-\sum _{i\in \mathcal{I}}{\lambda }_i{c}_i(x)-\frac{1}{2{\mu }_k}\sum _{i\in \mathcal{I}}({\lambda }_i-{\lambda }_i^k{)}^2\right\}.$$其中最后一项会对$\lambda$任何偏离前一次估计${\lambda }^k$的移动做出惩罚, 即, 它会使新的乘子$\lambda$与前一次估计${\lambda }^k$尽可能地相近. 因$\hat{F}$本身的定义就是一个关于$\lambda$的界约束二次规划, 因此我们有显式解$${\lambda }_i=\{\begin{array}{ll}0,& -c_i(x)+{\lambda }_i^k/{\mu }_k\le 0,\\ {\lambda }_i^k-{\mu }_k{c}_i(x),& otherwise.\end{array}$$注意到这可以用来更新乘子. 将此代入$\hat{F}$得到$$\hat{F}(x;{\lambda }^k,{\mu }_k)=f(x)+\sum _{i\in \mathcal{I}}\psi (c_i(x),{\lambda }_i^k;{\mu }_k),$$其中函数$\psi$具有三个标量参数:$$\psi (t,\sigma ;\mu )\stackrel{def}{=}\{\begin{array}{ll}-\sigma t+\frac{\mu }{2}{t}^2,& t-\sigma /\mu \le 0,\\ -\frac{1}{2\mu }{\sigma }^2,& otherwise.\end{array}$$因此, 对$x$极小化$\hat{F}(x;{\lambda }^k,{\mu }_k)$的$x_k$, 并依前面的显式表达更新Lagrange乘子. 比较框架3, 我们会发现$F$扮演了${\mathcal{L}}_A$的角色. 本节中介绍的格式就是等式约束增广Lagrange函数法向不等式约束情形的推广. 不过, 由于其性质还未得到验证, 此无约束化重构并不如界约束重构和约束线性化重构一般, 具有常用软件包的内嵌.

5. 不同算法的比较与软件介绍

在约束数目较小时, 人们通常使用二次罚函数法. 事实上, 有时仅需对一个较大的$\mu$做一次$Q(x;\mu )$的极小化. 若$\mu$选得不够好, 得到的解就可能不会很精确. 由于求解无约束优化的主要软件并不内嵌二次罚函数法, 因此很少有关于更新惩罚因子、调整容忍限以及选取初始点的讨论.

尽管二次罚函数更直观、简洁, 但第3节与第4节的增广Lagrange函数法通常更受青睐. 子问题的求解通常并不太难, 且乘子估计的引入避免了病态子问题的出现. 不过二次罚函数仍然作为许多算法正则化(regularization)的一种重要手段, 譬如SQP.

一般情形的$l_1$罚函数法由Fletcher在上世纪八十年代提出. 由于它与SQP算法的共同点, 它也被称为S$l_1$QP算法.最近, 软件包KNITRO加入了使用线性规划子问题的$l_1$罚函数法. 这两种方法我们将在下一章讨论.

近些年$l_1$罚函数受到了大量的关注. 它已被成功用于求解一些困难的问题, 例如带互补约束的数学规划(mathematical programs with complementarity constraints, MPCCs), 其中约束不满足标准的约束规范. 将这些问题约束吸收为罚项而不是线性化, 之后后再使用诸如SQP或内点法, 我们能够扩展这些其他方法的应用面. 软件包SNOPT在SQP算法中使用$l_1$罚函数法作为一种防护策略, 以防二次模型不可行或无界, 亦或有无界乘子.

增广Lagrange函数法由于其简洁性, 已受多年关注. MINOS和LANCELOT则是实施增广Lagrange函数法最好的软件包. 它们也适用于大型非线性规划问题. 一般层面上, MINOS的LCL和LANCELOT的BCL有共同的性质. 但它们在计算迭代步的子问题和求解子问题的技术上具有较大差异:

• MINOS使用既约空间法处理线性化约束, 并且使用(稠密)拟牛顿近似代替Lagrange函数的Hessian. 因此, MINOS对具有较小自由度的问题是最合适的.
• LANCELOT则更适用于约束较少的情形. 正如第4节所说, LANCELOT无需约束Jacobi矩阵$A$的分解, 这也增强了它在大型问题上的应用能力. 它同时也提供多种Hessian近似和预处理子的选择.
• PENNON软件包也基于增广Lagrange函数法, 它在处理半定矩阵约束时较有优势.

界约束重构和无约束重构Lagrange函数法的缺陷是, 它们增加了约束的平方, 让问题变得更复杂. 这时, 我们只有通过极小化增广Lagrange函数才能获得可行性的改善. 相反, LCL则在约束上计算牛顿类的迭代步, 稳步改善可行性. 因此, 不难预见, 在带线性约束的问题上MINOS比LANCELOT更具优势.

从第3节的增广Lagrange函数也可构造光滑精确罚函数(smooth exact penalty functions), 但这些就过于复杂了. 例如我们提过等式约束的Fletcher光滑精确罚函数:$${\varphi }_F(x;\mu )=f(x)-\lambda (x{)}^{T}c(x)+\frac{\mu }{2}\sum _{i\in \mathcal{E}}{c}_i(x{)}^2.$$其中Lagrange乘子估计$\lambda (x)$为最小二乘近似:$$\lambda (x)=[A(x)A(x{)}^T{]}^{-1}A(x)\nabla f(x)={\left(A^T\right)}^{†}\nabla f(x).$$这里的${\varphi }_F$光滑且精确. ${\varphi }_F$的缺点在于需要计算$\lambda (x)$. 当$A(x)$不满秩时$\lambda (x)$不是唯一的, 且当$A(x)$接近奇异时$\lambda$的估计可能会变得很差.

---

1. Han S P , Mangasarian O L . Exact penalty functions in nonlinear programming[J]. Mathematical Programming, 1979, 17(1):251-269. ↩︎

2. Han S P , Mangasarian O L . Exact penalty functions in nonlinear programming[J]. Mathematical Programming, 1979, 17(1):251-269. ↩︎