-
从1984年,Karmarkar[1]提出内点算法以来,内点算法的研究一直是一个热门的课题.随着研究的不断深入,研究者提出了许多不同类型的内点算法.从初始点的选取不同可分为:可行内点算法和不可行内点算法.不可行内点算法的优势是:不要求初始点严格可行,仅要求在非负象限的内部.由于这种类型的点容易获得,故在实际计算中显得尤为重要.
Lustig[2]首次提出了不可行内点算法,随后Masakazu[3]提出了原始-对偶不可行内点算法.这些早期的著作均没有给出不可行内点算法的多项式复杂度,直到Zhang[4]给出了第一个具有多项式复杂度O(n2 L)的不可行算法,然而,这个算法的复杂度较高.为了降低复杂度,众多学者做了大量研究工作[5-10].目前,不可行内点算法的最低复杂度为O(nL),这一结果仅对窄邻域不可行内点算法成立,并且窄邻域不可行内点算法的一般复杂度为O(n1.5L).对宽邻域而言一般算法的复杂度为O(n2L),目前最低复杂度为O(n1.5L).
本文将给出一个复杂度为O(n1.5L)的宽邻域不可行内点算法,该算法迭代在Ai[11]等提出的宽邻域内.最后,通过数值试验说明该算法是有效的.
本文中:e表示分量全为1的列向量;‖·‖表示向量的2 -范数.对于向量x,s∈
$ {{\mathbb{R}}^\mathit{n}} $ ,xs∈$ {{\mathbb{R}}^\mathit{n}} $ 表示对应分量的乘,(xs)min表示xs∈$ {{\mathbb{R}}^\mathit{n}} $ 的最小分量;让$ \mathit{\boldsymbol{p}} = {\mathit{\boldsymbol{x}}^{ - \frac{1}{2}}}{\mathit{\boldsymbol{s}}^{\frac{1}{2}}}, \mathit{\boldsymbol{\tilde r}} = {\left( {\mathit{\boldsymbol{xs}}} \right)^{\frac{1}{2}}}\mathit{\boldsymbol{r}} $ ;另外,对于h∈$ {{\mathbb{R}}^\mathit{n}} $ ,我们记:h+=max{0,h}和h-=h-h+.
全文HTML
-
本文将研究标准形式的原始-对偶线性规划:
其中:A∈
$ {{\mathbb{R}}^{m\times n}} $ ;c,x,s∈$ {{\mathbb{R}}^\mathit{n}} $ ;b,y∈$ {{\mathbb{R}}^\mathit{m}} $ .为获得(P)和(D)的近似解,一般将原始-对偶线性规划转化为扰动系统:
如果(P)和(D)的严格可行集
并且A行满秩,即rank(A)=m,则(1) 式存在唯一解.原始-对偶内点算法的基本思想是用Newton方法求解(1) 式,逐渐减小μ,并使得迭代点列包含在中心路径的某个邻域内,最终得到满足允许精度的最优解.本文假设F0≠
$ \phi $ 且rank(A)=m.下面给出求解迭代方向(Δx,Δy,Δs)的方程如下:
其中:
设α表示沿方向(Δx,Δy,Δs)的步长,则新的迭代为:
直接计算,可得下列常用结果:
本文将限制中心路径在宽邻域
$ \mathscr{N} $ (τ,β)内,其中$ \beta \le \frac{1}{2}, \tau \le \frac{1}{4} $ 是邻域参数并且另外,由(5) 式易知:
最后,本文通过下列方法计算最大步长
$ \tilde{\alpha } $ :
-
在前面准备工作的基础上,给出不可行内点算法的一般框架,如下:
初始化:
$ \beta \le \frac{1}{2}, \tau \le \frac{1}{4}, \left( {{\mathit{\boldsymbol{x}}}^{0}}, {{\mathit{\boldsymbol{y}}}^{0}}, {{\mathit{\boldsymbol{s}}}^{0}} \right)\in \mathscr{N}\left( \tau, \beta \right), {{\mu }^{0}}=\frac{1}{n}{{\left( {{\mathit{\boldsymbol{x}}}^{0}} \right)}^{\rm{T}}}{{\mathit{\boldsymbol{s}}}^{0}} $ .置k:=0.步骤1 如果(xk)Tsk≤ε(x0)Ts0,终止.
步骤2 通过解方程组(2) 获得方向(Δx,Δy,Δs),并通过(7) 式计算步长
$ {{\tilde{\alpha }}^{k}} $ .令步骤3 计算
$ {\mu ^{k + 1}} = \frac{1}{n}{\left( {{\mathit{\boldsymbol{x}}^{k + 1}}} \right)^{\rm{T}}}{\mathit{\boldsymbol{s}}^{k + 1}} $ ,并置k:=k+1,转到步骤1.给出下面关于本文算法的两点说明对分析算法的复杂度至关重要.
注1 若{(xk,yk,sk)}是上述算法产生的迭代序列,则对于k≥0,有
其中
从注1易知:
这蕴含着νk可以控制在点(xk,yk,sk)处的不可行性.
另外,由式(6) 可得:
这能够保证(xk)Tsk小于收敛精度时,不可行也可以小于相同的精度,等价于:当
$ \frac{{{{\left( {{\mathit{\boldsymbol{x}}^k}} \right)}^{\rm{T}}}{\mathit{\boldsymbol{s}}^k}}}{{{{\left( {{\mathit{\boldsymbol{x}}^0}} \right)}^{\rm{T}}}{\mathit{\boldsymbol{s}}^0}}} \le \varepsilon $ 时,必有注2 本文将选取与文献[12]中相同的初始点.
-
引理1[11] 若(x,y,s)∈
$ \mathscr{N} $ (τ,β),则引理2[13] 若(x,y,s)∈
$ \mathscr{N} $ (τ,β),则其中
由于本文使用了与文献[12]相同的初始点,故使用同样的证明方法可得下面引理.
引理3[12] 若
$ \left( {\mathit{\boldsymbol{x}}, \mathit{\boldsymbol{y}}, \mathit{\boldsymbol{s}}} \right) \in \mathscr{N}\left( {\tau, \beta } \right), \mathit{\boldsymbol{p}} = {\left( {{\mathit{\boldsymbol{x}}^{ - 1}}\mathit{\boldsymbol{s}}} \right)^{\frac{1}{2}}} $ ,则引理4 若
$ \left( {\mathit{\boldsymbol{x}}, \mathit{\boldsymbol{y}}, \mathit{\boldsymbol{s}}} \right) \in \mathscr{N}\left( {\tau, \beta } \right), \mathit{\boldsymbol{p}} = {\left( {{\mathit{\boldsymbol{x}}^{ - 1}}\mathit{\boldsymbol{s}}} \right)^{\frac{1}{2}}} $ ,则证 在方程(2) 的第三个等式两边同时乘以
$ {\left( {\mathit{\boldsymbol{xs}}} \right)^{ - \frac{1}{2}}} $ ,可得:进一步,可得:
其中第二个不等式使用了引理2.
另外,使用引理1和引理4容易获得下面的引理.
引理5 若(x,y,s)∈
$ \mathscr{N} $ (τ,β),则μ(α)在区间[0, 1]上是单调递减的.证 对(4) 式两边关于α求导可得:
这蕴含着μ(α)在区间[0, 1]上是单调递减的.
引理6 若(x,y,s)∈
$ \mathscr{N} $ (τ,β),则其中
证 使用(4) 式及引理3,可得:
另外,使用(4) 式及引理1、引理3可得:
证毕.
下面,将寻找满足条件(6) 的αf的下界.
引理7 若(x,y,s)∈
$ \mathscr{N} $ (τ,β)满足条件(6),则证 要证明该引理,只需要证明在区间
$ \left[{0, {{\tilde \alpha }^0}} \right] $ 上,x(α)Ts(α)≥(1-α)xTs成立.为此,首先证明其中最后一个不等式使用了引理3.
使用上面的不等式,可得:
证毕.
为了找到最大迭代步长
$ \tilde \alpha $ 的下界,首先需要给出下列引理.引理8[13] 若(x,y,s)∈
$ \mathscr{N} $ (τ,β),则:如果
$ \alpha \ge \frac{1}{{\sqrt n }} $ ,则如果
$ \alpha < \frac{1}{{\sqrt n }} $ ,则引理9 若(x,y,s)∈
$ \mathscr{N} $ (τ,β),$ {\tilde \alpha } $ 满足条件(7),则$ \tilde \alpha \ge {\tilde \alpha ^0} $ .证 为了证明
$ \left\| {\;{{\left( {\tau \mu \left( \alpha \right)\mathit{\boldsymbol{e}} - \mathit{\boldsymbol{x}}\left( \alpha \right)\mathit{\boldsymbol{s}}\left( \alpha \right)} \right)}^ + }\;} \right\|\; \le \beta \tau \mu \left( \alpha \right) $ ,首先,需要对不等式的左边进行放大,可得:其中第一个不等式使用了文献[11]的注解3.1,第二、第三个不等式分别使用了引理8、引理3.
构造一个函数f(α),将该引理的证明转化为证明f(α)≤0即可,其中
使用引理6以及
$ \beta \le \frac{1}{2}, \tau \le \frac{1}{4}, \omega \ge 6 $ ,可得:引理证毕.
下面给出算法的多项式复杂度.
定理1 让
$ \left( {{\mathit{\boldsymbol{x}}^0}, {\mathit{\boldsymbol{y}}^0}, {\mathit{\boldsymbol{s}}^0}} \right) \in {\mathscr{N}}\left( {\tau, \beta } \right), \tau \le \frac{1}{4}, \beta \le \frac{1}{2} $ ,则本文的算法至多需要O(n1.5L)次迭代停止,其中证 由算法的终止条件步骤1及引理6,可得:
其中第三个不等式使用了引理5,
又由于nμ0=(x0)Ts0,log(1-α)≤-α,容易证明
另外,当
$ k \ge \frac{{{n^{1.5}}}}{\eta }\log \left( {\frac{{{{\left( {{\mathit{\boldsymbol{x}}^0}} \right)}^{\rm{T}}}{\mathit{\boldsymbol{s}}^0}}}{\varepsilon }} \right) $ 时,使用(8) 式可得:所以,算法至多需要O(n1.5L)次迭代停止.
-
用本文算法(算法1) 和文献[6](算法2) 对文献[14]中给出的线性规划进行测试.用MATLAB R2011b编程并在Intel Core i3 PC (3.00GHz),2GB内存下运行.最终的测试结果见表 1.从表 1中容易发现:本文提出的算法在迭代次数和时间上均要优于文献[6]中提出的算法,尤其是对于大规模问题,这种优势更明显.