最优控制理论 六 拉格朗日乘子法和KKT条件

拉格朗日乘子法和KKT条件

1. 等式约束最优化2. 不等式约束最优化2.1 1个不等式约束2.2 KKT条件2.3 二维不等式约束图解3. MATLAB不等式约束优化总结4. 参考文献

最优控制是建立在最优化基础上的，它所处理的是无穷维路径函数的泛函极值问题，而后者是处理的对有限个参数所构成的函数的优化问题。也就是说，最优化处理的是min⁡xL(x),x∈Rn\min_{\mathbf x}L(\mathbf x),\mathbf x\in\Reals^nminx​L(x),x∈Rn；而最优控制处理的是min⁡x(t)L(x(t)),x(t)∈Rn\min_{\mathbf x(t)}L(\mathbf x(t)),\mathbf x(t)\in\Reals^nminx(t)​L(x(t)),x(t)∈Rn的问题，其中ttt为独立变量，x\mathbf xx都是被优化的变量。

对于最优控制问题而言，在处理等式约束条件时，采用的协态变量实际上就是时变的Lagrange乘子。具体而言，在变分法章节讲述终端约束ψ(x(tf),tf)=0\psi(x(t_f),t_f)=0ψ(x(tf​),tf​)=0时，直接引入了Lagrange乘数μ∈Rm\mu\in\Reals^mμ∈Rm；同样在处理路径等式约束f(x(t),x˙(t),t)=0f(x(t),\dot x(t),t)=0f(x(t),x˙(t),t)=0时，在Hamilton函数中引入了Lagrange乘子λ(t)∈Rn\lambda(t)\in\Reals^nλ(t)∈Rn。事实上，在最优化问题min⁡xL(x),s.t.g(x)=0\min_{\mathbf x}L(\mathbf x),s.t. g(\mathbf x)=0minx​L(x),s.t.g(x)=0 时，这种方法也是类似的。

本篇博客从等式约束最优化开始，进一步解释不等式约束最优化问题的必要条件，即KKT条件。

1. 等式约束最优化

下面举一个例子来说明拉格朗日乘子法在最优化问题中的应用：

例1：等式约束凸二次规划

min⁡x∈Rn12xTPx+qTx,s.t.Ax=b\min_{x\in\Reals^n}\frac 1 2x^\mathrm TPx+q^\mathrm Tx,\quad s.t.Ax=b x∈Rnmin​21​xTPx+qTx,s.t.Ax=b

其中b∈Rm,A∈Rn∗m,rank(A)=m,m≤nb\in\Reals^m,A\in\Reals^{n*m},rank(A)=m,m\leq nb∈Rm,A∈Rn∗m,rank(A)=m,m≤n。同样引入Lagrange乘子λ∈Rm\lambda\in\Reals^mλ∈Rm，构建Hamilton函数

H(x,λ)=12xTPx+qTx+λT(Ax−b)H(x,\lambda)=\frac 1 2x^\mathrm TPx+q^\mathrm Tx+\lambda^\mathrm T(Ax-b) H(x,λ)=21​xTPx+qTx+λT(Ax−b)

于是，等式约束优化问题转化为对无约束Hamilton函数的优化问题。最优解x∗x^*x∗和乘数λ∗\lambda^*λ∗服从一阶必要条件

∂H∂x∗=Px∗+q+Aλ∗=0∂H∂λ∗=Ax∗−b=0\frac{\partial H}{\partial x^*}=Px^*+q+A\lambda^*=0\\ \frac{\partial H}{\partial \lambda^*}=Ax^*-b=0 ∂x∗∂H​=Px∗+q+Aλ∗=0∂λ∗∂H​=Ax∗−b=0

该问题构成m+nm+nm+n维的线性方程组，若矩阵不奇异，则可以求出最优解。该问题的数值求解方法可见 zealscott-等式约束优化，不会的可以直接调用MATLAB函数fmincon\tt fminconfmincon.

利用**Lagrange乘子法，可以把等式约束问题转化为无约束问题，**这个无约束问题求解Hamilton函数的极值。下面我们再来了解一下Lagrange乘子的物理意义。问题 知乎-如何理解拉格朗日乘子法？从各个角度都进行了回答，我在这里简单摘取别人的话。

从平面几何角度。极值点处，函数等势面f(x)=cf(x)=\text cf(x)=c和等式约束g(x)=0g(x)=0g(x)=0的切线相切，也就是说两个梯度∇f(x),∇g(x)\nabla f(x),\nabla g(x)∇f(x),∇g(x)共线，而这个梯度前面的系数是Lagrange乘数λ\lambdaλ。梯度共线，意味着每个方向xix_ixi​上都有偏导数成比例，即

∂f∂xi=−λi∂g∂xi\frac{\partial f}{\partial x_i}=-\lambda_i\frac{\partial g}{\partial x_i} ∂xi​∂f​=−λi​∂xi​∂g​

那么写成向量形式就是∇f(x)+λT∇g(x)=0\nabla f(x)+\lambda^\mathrm T\nabla g(x)=0∇f(x)+λT∇g(x)=0

从代数角度。我们所定义的多元标量函数H(x,λ)=f(x)+λTg(x)H(x,\lambda)=f(x)+\lambda^\mathrm Tg(x)H(x,λ)=f(x)+λTg(x)，是一个无约束优化问题，它在极值点的必要条件就是，对自变量λ,x\lambda,xλ,x的梯度都为0。即

∂H∂x=∇f+λT∇g=0∂H∂λ=g(x)=0\frac{\partial H}{\partial x}=\nabla f+\lambda^\mathrm T\nabla g=0\\ \frac{\partial H}{\partial \lambda}=g(x)=0 ∂x∂H​=∇f+λT∇g=0∂λ∂H​=g(x)=0

第一个式子是必要条件，第二个式子是等式约束，这两个同时满足才行。

从函数构造方面来说，如果等式约束自然满足g(x)=0g(x)=0g(x)=0，那显然哈密尔顿函数的第二项就没有了，即H(x,λ)=f(x)H(x,\lambda)=f(x)H(x,λ)=f(x)。从这个意义上来说，所构造Hamilton函数和原问题是等价的。

2. 不等式约束最优化

上面讲述了等式约束最优化的处理方法，采用Lagrange乘子法把约束问题转化为无约束问题。下面我们来讨论不等式约束最优化问题，不同之处在于不等式约束是否active是需要分类讨论的。

2.1 1个不等式约束

为简单起见，首先考虑一个不等式的情况，例如g(x)=−x+c≤0g(x)=-x+c\leq0g(x)=−x+c≤0。

min⁡xf(x),s.t.g(x)≤0\min_xf(x),\quad s.t.g(x)\leq 0 xmin​f(x),s.t.g(x)≤0

不等式约束的处理方法要从不等式约束有没有起作用（active）谈起。阴影部分为与不等式冲突的区域，即g(x)>0g(x)>0g(x)>0，下图反映了不等式与极小值之间的两种情况

图1，最小值本来就在g(x∗)<0g(x^*)<0g(x∗)<0的部分取得；图2，由于受到不等式的约束，函数极小值只能在边界处取得，此时，g(x∗)=0g(x^*)=0g(x∗)=0。如果引入松弛变量μ∈R\mu\in\Realsμ∈R，把不等式约束暂且用等式约束的方法来处理，构建如下形式的Hamilton函数，H(x,μ)=f(x)+μg(x)H(x,\mu)=f(x)+\mu g(x)H(x,μ)=f(x)+μg(x)。函数取得极值的一阶必要条件为

∂H∂x=∇f(x∗)+μ∇g(x∗)=0g(x∗)≤0μ{=0,g(x∗)<0(本来就满足,约束条件不影响结果)>0,g(x∗)=0(限制了最小值的取得)\begin{aligned} &\frac{\partial H}{\partial x}=\nabla f(x^*)+\mu \nabla g(x^*)=0\\ &g(x^*)\leq 0\\ &\mu\begin{cases} =0, &g(x^*)<0&\text{(本来就满足,约束条件不影响结果)}\\ \gt0, &g(x^*)=0&\text{(限制了最小值的取得)} \end{cases} \end{aligned} ​∂x∂H​=∇f(x∗)+μ∇g(x∗)=0g(x∗)≤0μ{=0,>0,​g(x∗)<0g(x∗)=0​(本来就满足,约束条件不影响结果)(限制了最小值的取得)​​

上式∇(⋅)=∂(⋅)∂x\nabla (\cdot)=\frac{\partial(\cdot)}{\partial x}∇(⋅)=∂x∂(⋅)​，把两种情况写作一个条件，即

{μ≥0μ⋅g(x∗)=0\left\{\begin{matrix} \mu\geq 0\\ \mu\cdot g(x^*)=0 \end{matrix}\right. {μ≥0μ⋅g(x∗)=0​

下面解释为什么松弛变量μ≥0\mu\geq 0μ≥0而不能小于0. 我们引入松弛变量μ\muμ的目的是，使新的问题和原始问题等价，只有这样才能保证等价。大于号和等于号对应了不等式约束对优化问题的限制行为，这个由图1和图2就可以明显看出来。当μ=0\mu=0μ=0时，不等式约束是可有可无的，没有它求出的极小值和原来的一样。当μ>0\mu>0μ>0时，和约束冲突的解∀x,g(x)>0\forall x,g(x)>0∀x,g(x)>0必然使得H(x,μ)=f(x)+μg(x)>f(x)H(x,\mu)=f(x)+\mu g(x)>f(x)H(x,μ)=f(x)+μg(x)>f(x)，这也就避免了它取这个点。从一阶必要条件来看，如果极小值就恰好在不等式的边界取得，此时，∇f(x)>0(→),∇g(x)<0(←)\nabla f(x)>0(\rightarrow),\nabla g(x)<0(\leftarrow)∇f(x)>0(→),∇g(x)<0(←)，要使两者达到平衡的系数是μ>0\mu>0μ>0，这样∇f(x∗)+μ∇g(x∗)=0\nabla f(x^*)+\mu \nabla g(x^*)=0∇f(x∗)+μ∇g(x∗)=0.

2.2 KKT条件

2.1节函数极小值的必要条件即为KKT条件，Karush–Kuhn–Tucker conditions。参考doubleslow- KKT条件详解这篇文章，以下直接给出KKT条件。设不等式+等式约束的非线性最优化问题如下

min⁡xf(x),s.t.h(x)=0g(x)≤0\min_\mathbf x f(\mathbf x),\quad s.t.\begin{aligned} \mathbf h(\mathbf x)=0\\ \mathbf g(\mathbf x)\leq0 \end{aligned} xmin​f(x),s.t.h(x)=0g(x)≤0​

其中x∈Rn\mathbf x\in\Reals^nx∈Rn为多元向量，f(x)f(x)f(x)为标量函数，h∈Rp,g∈Rq\mathbf h\in\Reals^p,\mathbf g\in\Reals^qh∈Rp,g∈Rq为多元向量函数，对该问题引入Lagrange乘子λ∈Rp\lambda\in\Reals^pλ∈Rp，松弛变量μ∈Rq\mu\in\Reals^qμ∈Rq. 构建Hamilton函数

H(x,λ,μ)=f(x)+λTh(x)+μTg(x)H(x,\lambda,\mu)=f(x)+\lambda^\mathrm T h(x)+\mu^\mathrm Tg(x) H(x,λ,μ)=f(x)+λTh(x)+μTg(x)

该问题的KKT条件为

{∇H=0(一阶必要条件)h(x∗)=0,g(x)≤0(原始约束)λi≠0(等式约束始终起作用)μjgj(x∗)=0(互补松弛性)μj≥0(松弛变量非负)\left\{\begin{matrix} \nabla H=0& (一阶必要条件)\\ \mathbf h(\mathbf x^*)=0, \mathbf g(\mathbf x)\leq0 &(原始约束)\\ \lambda_i\neq0&(等式约束始终起作用)\\ \mu_jg_j(\mathbf x^*)=0&(互补松弛性)\\ \mu_j\geq0&(松弛变量非负) \end{matrix}\right. ⎩⎪⎪⎪⎪⎨⎪⎪⎪⎪⎧​∇H=0h(x∗)=0,g(x)≤0λi​​=0μj​gj​(x∗)=0μj​≥0​(一阶必要条件)(原始约束)(等式约束始终起作用)(互补松弛性)(松弛变量非负)​

上式∇(⋅)=∂(⋅)∂x∈Rn，i=1,2,⋯,p,j=1,2,⋯,q\nabla (\cdot)=\frac{\partial(\cdot)}{\partial x}\in\Reals^n，i=1,2,\cdots,p, j=1,2,\cdots,q∇(⋅)=∂x∂(⋅)​∈Rn，i=1,2,⋯,p,j=1,2,⋯,q.

2.3 二维不等式约束图解

下面是一个不等式约束的非线性规划问题，

min⁡x,yf(x,y)=xexp⁡[−(x2+y2)]+(x2+y2)/20s.t.g(x,y)=xy/2+(x+2)2+(y−2)2/2−c≤0\min_{x,y}f(x,y)=x \exp[-(x^2+y^2)]+(x^2+y^2)/20\\ s.t.g(x,y)=xy/2+(x+2)^2+(y-2)^2/2-c\leq0 x,ymin​f(x,y)=xexp[−(x2+y2)]+(x2+y2)/20s.t.g(x,y)=xy/2+(x+2)2+(y−2)2/2−c≤0

同样用梯度法来解释，有下面的二维图形，红色弧线内部是不等式满足的区域，即g(x)<0g(x)<0g(x)<0。图3的极小值点是无约束问题本来的极小值点。在这一点x∗x^*x∗，不等式约束不起作用，要使得H(x,y,μ1,μ2)=f(x,y)H(x,y,\mu_1,\mu_2)=f(x,y)H(x,y,μ1​,μ2​)=f(x,y)，所以μj=0\mu_j=0μj​=0。图4中，极小值是在不等式边界上取得的。在这一点x∗x^*x∗，梯度沿x,yx,yx,y方向互相成比例，即∇xf(x∗)+μ1∇xg(x∗)=0,∇yf(x∗)+μ2∇yg(x∗)=0\nabla_x f(x^*)+\mu_1 \nabla_x g(x^*)=0,\nabla_y f(x^*)+\mu_2 \nabla_y g(x^*)=0∇x​f(x∗)+μ1​∇x​g(x∗)=0,∇y​f(x∗)+μ2​∇y​g(x∗)=0 ，由于梯度是反向的，所以μj>0\mu_j>0μj​>0.

实际上和下面这篇文章的表达方式一样，知乎-图解KKT条件和拉格朗日乘子法 ，绘图解释了二维函数+线性不等式约束KKT条件。

3. MATLAB不等式约束优化

很多形式规范的约束最优化问题都可以直接调用MATLAB中现有的优化函数，下面是各类情形常用来的求解算法：

需要注意的是，并不是所有优化问题都可以用这些算法，（比如一些有很多局部最优解的问题就只能用智能优化算法），以上表格中列举的算法求解对应形式的问题效率非常高，而采用智能算法的效率就不一定更高。

总结

一言以蔽之，等式约束的最优性条件就是说被优化函数的梯度，等于多个等式约束函数梯度的线性组合，由于等式约束是互相线性无关的，所以等式约束的函数梯度的线性组合系数就是拉格朗日乘子。

不等式约束的自由训条件是说被优化函数的梯度，等于多个等式约束以及起作用的(active)不等式约束的梯度的线性组合，然后这个组合系数的就是每个等式的线性

4. 参考文献

[1] Bryson A E , Ho Y C ,Applied optimal control : optimization, estimation, and control Sec 1.1~Sec 1.7 [J]. IEEE Transactions on Systems Man & Cybernetics, 1975

[2]. 见文中的链接。

[3]. 桂。cnblogs 最优化：拉格朗日乘子法