乘子法

乘子法时针对外部罚函数法的改进方法，由于外部罚函数法随着罚因子的增大，增广目标函数的Hesse矩阵条件数会逐渐增大，从而导致在实际计算中，数值计算的稳定性也会变得越来越差，难以精确求解。乘子法改由在约束问题的Lagrange函数中加入相应的惩罚，使得在求解无约束问题时，罚因子不必趋于无穷大就能求到约束问题的最优解，而且数值计算的稳定性也能得到很好的保证。

等式约束

等式约束的H乘子法

H乘子法全称为Hesternes乘子法。
首先考虑如下的等式约束问题(用向量形式进行表示)

$$\left.\begin{array}{rl}min&f(x)\\s.t.&\boldsymbol h(x)=0\end{array}\right\}\tag1$$ 其中 $f,\boldsymbol h$时二次连续可微函数，设 $D=\{x|\boldsymbol h(x)=0\}$，即 $D$表示容许集。
那么问题(1)的Lagrange函数为 $$L(x,\lambda)=f(x)-\boldsymbol \lambda^T\boldsymbol h(x)\tag2$$ 假设 $x^*$时问题(1)的极小点， $\lambda^*$时相应的Lagrange乘子，那么必有 $$\nabla_x L(x^*,\lambda^*)=\nabla f(x^*)-\nabla \boldsymbol h(x^*)\lambda^*=0\tag3$$ $$\nabla_{\lambda}L(x^*,\lambda^*)=h(x^*)=0\tag4$$ 对于 $\forall x\in D$都有 $\boldsymbol h(x)=0$，因此 $$f(x^*)=L(x^*,\lambda^*)\le L(x,\lambda^*)=f(x)\tag5$$ 所以问题(1)可以等价为 $$\left.\begin{array}{rl}min&L(x,\lambda^*)\\ s.t.&\boldsymbol h(x)=0\end{array}\right\}\tag6$$ 也就是说我们可以通过公式(6)来求解问题(1)。公式(6)的增广目标函数为 $$F(x,\lambda^*,\mu)=L(x,\lambda^*)+\mu\boldsymbol h(x)^T\boldsymbol h(x)\tag7$$ 但是我们知道公式(7)中的 $\lambda^*$其实时未知的。而乘子法的基本思想就是在求解 $x^*$的同时，通过迭代的方法，同时求解出 $\lambda^*$。

首先需要明确的是，使用乘子法求解问题(1)时的增广目标函数是

$$F(x,\lambda,\mu)=f(x)-\boldsymbol \lambda^T\boldsymbol h(x)+\mu\boldsymbol h(x)^T\boldsymbol h(x)\tag8$$

下面我们会直接给出一些引理和定理，不会对引理和定理进行证明。
引理1 设$U$是$n\times n$对称矩阵，$H$是$n\times l$矩阵，如果满足$H^Tx=0$的所有非零向量$x$都有$x^TUx\gt0$，那么一定存在非负数$\mu^*$，当$\mu\ge\mu^*$时，使得矩阵$U+\mu HH^T$在$R^n$上都正定。
定理2 在约束问题(1)中，如果最优解$x^*\in R^n$和Lagrange乘子$\lambda^*\in R^l$满足二阶充分条件(见最优性条件的定理2)，那么就存在非负数$\mu^*$，使得当$\mu\ge\mu^*$时，$x^*$也是$F(x,\lambda^*,\mu)$的严格局部极小点。
该定理表明，对公式(6)使用外部罚函数法的时候，只要罚因子不小于某一个潜在的非负数，就能保证$F(x,\lambda^*,\mu)$的(局部)极小点恰是(1)的(局部)极小点。
现在求解以(8)为目标函数的无约束极小化问题。

$$min\quad F(x,\lambda,\mu)\tag9$$
定理3 设对于给定的参数 $\lambda$和 $\mu$， $x^*$是公式(9)的最优解，那么 $x^*$是问题(1)的最优解，且 $\lambda$为相应的Lagrange乘子向量 $\lambda^*$的充要条件是 $\boldsymbol h(x^*)=0$，即 $x^*\in D$。
定理2结合定理3说明:如果恰好选取参数 $\lambda$为 $\lambda^*$，选择 $\mu\ge\mu^*$，然后求解公式(9)而得出的最优解 $x^*$恰好是问题(1)的容许点，那么 $x^*$就是问题(1)的最优解。

算法:等式约束问题的H乘子法
已知：目标函数$f(X)$，等式约束函数$h_1(x),h_2(x),...,h_l(x)$；终止限$\epsilon$；常数$\theta\in(0,1)$
(1)选定初始点$x_0$，乘子向量的初始估计$\lambda_1$，初始罚因子$\mu$，罚因子的放大系数$c(\gt1)$，置$k=1$
(2)以$x_{k-1}$为初始点，求解无约束问题$\min\limits_x\quad F(x,\lambda_k,\mu)$
(3)若$||\boldsymbol h(x_k)||\le\epsilon$，则打印$x_k$，停止；否则，转(4)
(4)若$||\boldsymbol h(x_k)||/||\boldsymbol h(x_{k-1})||\le\theta$，则转(5)；否则，置$\mu=c\mu$，然后转(5)
(5)计算$\lambda_{k+1}=\lambda_k-2\mu h(x_k)$；置$k=k+1$，转(2)。

通过以上算法我们知道，$\lambda^*$是通过迭代求出的，迭代公式为

$$\lambda_{k+1}=\lambda_k-2\mu h(x_k)$$ ，此迭代公式的推导如下。
对于给定的 $\lambda=\lambda_k$，设 $x_k$是公式(9)，即 $$\min\limits_xF(x,\lambda_k,\mu)\tag{10}$$ 的极小点，那么应有 $$\nabla_xf(x_k,\lambda_k,\mu)=\nabla f(x_k)-\nabla h(x_k)[\lambda_k-2\mu h(x_k)]=0\tag{11}$$ 问题(1)的最优点 $x^*$处，应有 $$\nabla_xL(x^*,\lambda^*)=\nabla f(x^*)-\nabla h(x^*)\lambda^*=0\tag{12}$$ 若设 $\nabla h(x^*)$是满秩矩阵，则使用最小二乘法，由上式可以求出 $$\lambda^*=[\nabla h(x^*)^T\nabla h(x^*)]^{-1}\nabla h(x^*)^T\nabla f(x^*)$$ 且若 $x_k\rightarrow x^*$，则当 $k$充分大时，$\nabla h(x_k)$也将是满秩矩阵，于是，再次使用最小二乘法，由公式(11)又可以得到 $$\lambda_k-a\mu h(x_k)=[\nabla h(x_k)^T\nabla h(x_K)]^{-1}\nabla h(x_k)^T\nabla f(x_k)$$ 比较以上两个式子克制，当 $x_k\rightarrow x^*$时， $\lambda_k-2\mu h(x_k)\rightarrow\lambda^*$，所以有理由采用迭代公式 $$\lambda_{k+1}=\lambda_k-2\mu h(x_k)$$

等式约束的P乘子法

P乘子法全称为Powell乘子法，与H乘子法的差别时它对每一个约束使用了不同的罚因子。P乘子法的增广目标函数为

$$F(x,\lambda,\mu)=f(x)-\sum\limits_{i=1}^l\lambda_ih_i(x)+\sum\limits_{i=1}^l\mu_ih_i^2(x)\tag{13}$$

算法：等式约束的P乘子法
已知：目标函数$f(x)$，等式约束函数$h_1(x),h_2(x),...,h_l(x)$；终止限$\epsilon$；常熟$\theta\in(0,1)$
(1)选定初始点$x_0$，乘子向量的初始估计$\lambda_1$，初始罚因子向量$\mu_1$，罚因子的放大系统$c$；置$k=1$
(2)以$x_{k-1}$为初始点，求解以公式(13)为目标函数的极小化问题

$$min F(x,\lambda_k,\mu_k)$$ 设其极小点为 $x_k$
(3)若 $||h(x_k)||\lt\epsilon$，则打印 $x_k$，停止；否则转(4)
(4)对 $j=1,2,...,l$，一次判别 $|\frac{h_i(x_k)}{h_j(x_{k-1})}|\lt\theta$是否成立：若成立，则置 $\mu_j^{(k+1)}=\mu_j^{(k)}$；否则，置 $\mu_j^{(k+1)}=c\mu_j^{(k)}$，直至判别结束，然后转(5)
(5)对 $j=1,2,...,l$一次计算 $\lambda_j^{(k+1)}=\lambda_j^{(k)}-2\mu_j^{(k+1)}h_j(x_K)$；置 $k+k+1$，转(2)

计算实践证明，P乘子法比H乘子法效果好。

不等式约束

这里将介绍Rockafellar乘子法，他是把H乘子法巧妙的用到不等式约束问题上而得到的求解仅含不等式约束问题的一种方法。
考虑不等式约束问题

$$\left.\begin{array}{rl}min&f(x)\\s.t.&s_i(x)\ge0\quad i=1,2,...,m\end{array}\right\}\tag{14}$$ 引入辅助变量 $z_1,z_2,...,z_m$，把公式(14)转换为等式约束问题 $$\left.\begin{array}{rl}min&f(x)\\s.t.&s_i(x)-z_i^2=0\quad i=1,2,...,m\end{array}\right\}\tag{15}$$ 这时可以使用P乘子法求解(15)，其增广目标函数为 $$\tilde F(x,z,v,\mu)=f(x)-\sum\limits_{i=1}^mv_i[s_i(x)-z_i^2]+\mu\sum\limits_{i=1}^m[s_i(x)-z_i^2]^2\tag{16}$$ 求 $\tilde F$的极小点，由于变量 $x，z$独立，所以求 $\tilde F$的极小可分步进行，即 $$\min\limits_{x,z}\tilde F(x,z,v,\mu)=\min\limits_x\{\min\limits_z\tilde F(x,z,v,\mu)\}$$ 实际上， $\tilde F$关于 $z$的绩效可以用$x,v,\mu$显示表出，因此将(14)的增广目标函数定义为 $$\begin{array}{rl}F(x,v,\mu)&=\min\limits_z\tilde F(x,z,v,\mu)\\&=f(x)+\sum\limits_{i=1}^m\min\limits_{z_i}\{-v_i[s_i(x)-z_i^2]+\mu[s_i(x)-z_i^2]^2\}\end{array}\tag{17}$$ 用解析法求解(17)，令 $$\nabla_z\tilde F(x,z,v,\mu)=0$$ 即 $$v_iz_i-2\mu[s_i(x)-z_i^2]z_i=0\quad i-1,2,...,m$$ 经整理得 $$z_i[2\mu z_i^2-(z\mu s_i(x)-v_i)]=0\quad i=1,2,...,m\tag{18}$$ 由(18)，对 $i=1,2,...,m$有 $$z_i^2=\begin{cases}0&若2\mu s_i(x)-v_i\le 0\\s_i(x)-\frac{v_i}{2\mu}&若2\mu s_i(x)\gt0\end{cases}$$ 因此 $$s_i(x)-z_i^2=\begin{cases}s_x(x)&若s_i(x)\le\frac{v_i}{2\mu}\\\frac{v_i}{2\mu}&若s_i(x)\gt\frac{v_i}{2\mu}\end{cases}$$ 进而 $$\begin{array}{l}-v_i[s_i(x)-z_i^2]+\mu[s_i(x)-z_i^2]^2$\\=\begin{cases}-v_is_i(x)+\mu s_i^2(x)&当s_i(x)\le\frac{v_i}{2\mu}\\-\frac{v_i^2}{2\mu}+\mu[\frac{v_i}{2\mu}]^2&当s_i(x)\gt\frac{v_i}{2\mu}\end{cases}\\=\begin{cases}\mu[s_i(x)-\frac{v_i}{2\mu}]^2&当s_i(x)\le\frac{v_i}{2\mu}\\-\frac{v_i^2}{4\mu}&当s_i(x)\gt\frac{v_i}{2\mu_i}\end{cases}\\=\frac1{4\mu}\{[max(0,v_i-z\mu s_i(x))]^2-v_i^2\}\end{array}$$ 带回到(17)中，得到 $$F(x,v,\mu)=f(x)+\frac1{4\mu}\sum\limits_{i=1}^m\{[max\{0,v_i-2\mu s_i(x)\}]^2-v_i^2\}$$ 乘子迭代公式为 $$v_{k+1}=v_k$-2\mu[s(x_k)-z_k^2]=max\{0,v_k-2\mu s(x_k)\}$$ 即 $$v_i^{k+1}=\max\{0,v_i^{(k)}-2\mu s_i(x_k)\}\quad i=1,2,...,m$$ 终止准则为 $$\left(\sum\limits_{i=1}^m\left[\min\left\{s_i(x_k),\frac{v_i^{(k)}}{2\mu}\right\}\right]^{1/2}\right)\lt\epsilon$$

一般约束

对于一般约束问题

$$\begin{array}{rl}min&f(x)\\s.t.&s_i(x)\ge0\quad i=1,2,...,m\\&h_j(x)=0\quad j=1,2,...,l\end{array}$$
可以得到增广目标函数为 $$F(x,v,\lambda,\mu)=f(x)+\frac1{4\mu}\sum\limits_{i=1}^m\{[\max\{0,v_i-z\mu s_i(x)\}]^2\}-\sum\limits_{j=1}^l\lambda_jh_i(x)+\mu\sum\limits_{j=1}^l[h_j(x)]^2\tag{19}$$
乘子迭代公式为 $$\begin{cases}v_i^{(k+1)}=\max\{0,v_i^{(k)-2\mu s_i(x_k)}\}\quad i=1,2,...,m\\ \lambda_j^{(k+1)}=\lambda_j^{(k)}-2\mu h_j(x_k)\qquad j=1,2,...,l\end{cases}\tag{20}$$ 其中 $\lambda_j$代表的是第 $j$个等式约束所对应的Lagrange乘子，$v_i$代表的是第 $i$个不等式约束所对应的Lagrange乘子，显然，$v_i\ge0$

算法:一般约束问题的H乘子法
已知：目标函数$f(x)$，不等式约束函数$s_1(x),...,s_m(x)$，等式约束函数$h_1(x),...,h_l(x)$；终止限$\epsilon$；尝试$\theta\in(0,1)$，设中驰判别函数为

$$\phi_k=\left(\sum\limits_{i=1}^m\left[\min\left\{s_i(x_k),\frac{v_i^{(k)}}{2\mu}\right\}\right]^2+\sum\limits_{j=1}^l[h_j(x_k)]^2\right)^{1/2}$$
(1)选定初始点 $x_0$，乘子向量初始估计 $v_1$和 $\lambda_1$，初始罚因子 $\mu$，罚因子方法系数 $c(\gt1)$；计算 $\phi_0$；置 $k=1$
(2)以 $x_{k-1}$为初始点，求解 $$\min\limits_x F(x,v_k,\lambda_k,\mu)$$ 设 $x_k$为其极小点。
(3)计算 $\phi_k$，若 $\phi_k\lt\epsilon$，则输出 $x_k$，停止；否则，转(4)
(4)若 $\phi_k/phi_(k-1)\le\theta$，则转(5)；否则，置 $\mu=c\mu$，转(5)
(5)按公式(20)计算乘子 $v_i^{(k+1)},...,v_m^{(k+1)},\lambda_1^{(k+1)},...,\lambda_l^{(k+1)}$，置 $k=k+1$,转(2)