划:mincTxsubject to x≥0 首先令 c=(11),t(0)=4,x 的初始点为(1.3,1)
下面两图是第一次迭代,t=4 时的 tf0(x)+(x) 函数图象和平面的等高图,同时第二张图上还标出了迭代点的轨迹。
下面是t=8 时的 tf0(x)+(x) 函数图象和平面的等高图,同时第二张图上还标出了迭代点的轨迹。
对于本例,central path 就是 x1=x2 这条直线。从迭代点的轨迹可以看出,在第一次迭代后,迭代点一直在central path上移动。
其实对于Barrier Method,我自己之前有个疑问就是,既然我们知道在 x(t) 处的目标函数值与最优解的差肯定小于
mt,那为什么不直接给一个较大的 t 值通过较少的几次迭代就能算出最有解了呢。真对这一问题,原因可能是这样的,Barrier Method的计算量是由内外两层迭代组成的,外层对 t 进行更新:t=μt,内层用牛顿法求 tf0+ 的解。有仿真结果表明,随着 μ 值的加大,总的计算量(牛顿迭代次数)在呈现一定阶段的下降后并没有明显下降。而对与 μ 值的选取则是一个内层迭代次数和外层迭代次数的权衡。μ 值较大时,外层迭代次数少,内层迭代次数多;μ 值较小时情况则相反。一般来说 μ 的取值在 10到20 之间比较合适。 原始对偶内点法(Primal Dual Interior Point Method) 在 Primal Dual 的方法中,我们直接考虑一个标准的 LP 命题。mincTx,subject to Ax=b,x≥0(9) 它的对偶命题为:maxbTλ,subject to ATλ+s=c,s≥0(10) 上面两个命题的KKT条件为:
ATλ+sAxxisi(x,s)=c,=b,=0,i=1,...,n≥0.(11a)(11b)(11c)(11d) 为了后面的推导方便,将上面的KKT条件化为矩阵形式如下:F(x,λ,s)=ATλ+scAxbXSe(x,s)=0.≥0,(12a)(12b) 其中,X=diag(x1,x2,...,xn),S=diag(s1,s2,...,sn),e=(1,1,...1)T 如同一般的优化算法,这里需要定义一个量来检验当前的迭代点与最优点的差距。在Barrier Method中,使用 duality gap 的上界 mt 来检验的,在 Primal Dual Method 中,我们定义一个新的 duality measure 来进行某种衡量:
μ=1n∑i=1nxisi=xTsn 注意:这里的 μ 与 Barrier Method 中的 μ 意义不同,为了与各自书中的表达方式对应,分别选用了各自书中的变量记法。
要求解原始的优化命题,需要做的就是去求解 (12) 这样一个方程组,由于 F 阵第三行中 XS 一项的存在,因此这是一个非线性系统。要求解这样一个非线性系统,一种实用的方法就是牛顿法。(注意,这里说的牛顿法是一种求解非线性方程组的方法,与求解优化命题的牛顿法并不完全一样,但核心思路是一致的,都是在迭代点处进行二阶展开。)在当前点处求解一个前进方向,并通过迭代逼近非线性系统的解。J(x,λ,s)ΔxΔλΔs=F(x,λ,s) 其中 J 是F 阵的雅各比矩阵。我们将 F 阵中的前两行分别记为:rb=Axb,rc=ATλ+sc(13) 那么在每次迭代中需要求解的线性系统为:
0ASAT00I0XΔxΔλΔs=rcrbXSe 因此,在求解得到相应的前进方向后,需要做的就是求解确定一个步长 α 来进行如下的更新 (x,λ,s)+α(Δx,Δλ,Δs) 其中 α∈(0,1] 的选取要使得原始变量和对偶变量都满足相应的约束。
看起来这种方法似乎已经可以了,但通过这种被称为 affine scaling direction 的方向所得到的 α 往往很小。也就是说,很难在迭代中取得进展。一开始看到这个说法时,我也很难理解这句话的意思。所以在这里我们用一个图来说明,令 (9) 中的 c=(101),初始点为(0.7,2),这里不考虑等式约束,直
接采用affine scaling direction 的方向得到的迭代点的轨迹为
从图中可以看出,几乎在从第二次迭代开始,迭代点就一头撞上了约束。后面的迭代也只能贴着约束的边缘来走,因为要保迭代点不能违反约束,因此每次的步长 α 都只能取得很小。在这一过程中,一共进行了11次迭代才使得 duality measure μ<105 。
因此,实际上内点法采用的是一个不那么 aggressive 的牛顿方向,也就是控制迭代点使其徐徐想约束边界和最优点处靠近。具体的方法是,我们在用牛顿法求解非线性系统时,在每次迭代中并不要求直接实现 xisi=0,而是令其等于一个逐渐减少的值,具体来说就是 xisi=σμ,其中 μ 是当前的 duality measure,σ∈[0,1] 是用于控制下降速度的下降因子。也就是说,在每次迭代中要求解的方程组应该是
0ASAT00I0XΔxΔλΔs=rcrbXSe+σμe(14) 其中,σ 被称为 center parameter 。意在表示我们要通过调整 σ 使得迭代点的轨迹在 central path 附近。 Central Path
在内点法中,central path 是指满足下面一组方程式的迭代
点所组成的所组成的一条弧线
ATλ+sAxxisi(x,s)=c,=b,=τ,i=1,2...n.>0(15a)(15b)(15c)(15d) 这与 (11) 中的KKT 条件的区别就在于在第三个条件的等式右端的 0 被 τ>0 替代了。
也就是说,对偶内点法的基本思路就是依据 central path 计算出相应的方向,并控制步长 α 的选择使得迭代点位于 central path 的某一个邻域内。(有关central path 的邻域的内容,和具体的原始对偶内点法的算法在这里不作详细说明,有兴趣这可以参考相应书籍。)
关于步长 α 的选取,central parameter σ 的更新,以及更多的收敛性分析的内容,在这里不作展开。 举例
与上一张图对应,我们同样取 c=(101),初始点为(0.7,2),不考虑等式约束。采用对偶内点法的迭代点的轨迹为
可以看出,引入 central path 后,迭代点能在贴近约束边界后再次远离约束边界(也就是靠近 central path) ,从而保证下一次迭代能取得更大的进步。在这一过程中,一共进行了6次迭代就使得 duality measure μ<105 。 几个问题