劉新澤,劉紅衛(wèi),劉長河

(1.西安電子科技大學 數(shù)學系,西安 710071；2.臨滄高等師范?？茖W校 數(shù)理系,云南 臨滄 677000；3.河南科技大學 數(shù)學與統(tǒng)計學院,河南 洛陽 471003)

0 引 言

自從Karmarkar[1]提出線性規(guī)劃的第一個內(nèi)點算法后,內(nèi)點算法即成為運籌學領域的研究熱點之一.內(nèi)點算法不僅形式簡潔,而且實際執(zhí)行非常有效.實踐結果表明,預估-校正內(nèi)點算法是求解線性規(guī)劃、線性互補問題的有效方法[2-3].內(nèi)點算法中大鄰域算法實際計算效果較好,但理論復雜性相對較差;而小鄰域算法則相反.文獻[4-6]降低了大鄰域算法的迭代復雜度;文獻[7-8]利用預估-校正策略改進了算法的計算效果.但這些算法都是可行內(nèi)點算法,即需要嚴格初始可行點.而在實際問題中,可能不存在嚴格可行點,或者不容易找到一個嚴格初始可行點.因此,內(nèi)點算法的執(zhí)行大多數(shù)采用不可行內(nèi)點算法.由于不可行算法的迭代點不滿足等式約束,因此算法的收斂性和復雜性分析變得十分困難.本文基于文獻[5]提出的大鄰域算法,給出一種求解非單調(diào)線性互補問題(LCP)的不可行二階預估-校正算法,該算法具有O((1+κ)5/2nL)復雜度.數(shù)值實驗驗證了算法的有效性.

LCP的一般形式是: 求(x,s)∈n×n,使得

x≥0,s=Ax+q≥0,xs=0,

(1)

其中:A:n→n是線性算子;q∈n.如果A滿足: 存在常數(shù)κ≥0,使得對任意的u∈n,都有其中:I+={i:ui(Au)i>0};I-={i:ui(Au)i<0},則稱問題(1)為P*(κ)-LCP.顯然,P*(κ)-LCP是單調(diào)互補問題的推廣.本文記問題(1)的可行集為F={(x,s):x≥0,s≥0,s=Ax+q},解集為F*={(x*,s*)∈F: (x*)Ts*=0}.假設F*≠?.約定: ‖x‖(‖x‖1)表示向量x∈n的2-范數(shù)(1-范數(shù));e=(1,…,1)T;x≥0(x>0)表示xi≥0(xi>0),i=1,2…,n,x+=max{x,0},x-=min{x,0},X=diag(x)為向量x生成的對角矩陣;xs∶=Xs.

1 算 法

假設(x*,s*)是問題(1)的任一解,則存在正數(shù)ξ>0,使得ξ≥‖(x*,s*)‖.算法的初始點為(x0,s0)=ξ(e,e).記q0=s0-Ax0.考慮問題(1)的擾動系統(tǒng):

Ax-s+q=τ(q-q0),xs=μe,

(2)

其中:τ∈(0,1];μ=xTs/n.設(x(τ),s(τ))=(1-τ)(x*,s*)+τ(x0,s0),τ∈(0,1].易驗證x(τ)>0,s(τ)>0是問題(1)的嚴格可行解.因此,對任意的(μ,τ)∈(0,∞)×(0,1],式(2)有唯一解(x(μ,τ),s(μ,τ)),這些解構成中心路徑[9].當(μ,τ)→(0,0)時,中心路徑的極限即為問題(1)的一個最優(yōu)解.內(nèi)點法用牛頓法求解方程組(2),逐漸減小μ,并使得迭代點列包含于中心路徑的某個鄰域內(nèi),最終得到滿足允許精度的最優(yōu)解.本文算法使用如下大鄰域[5]:

N(β,γ)={(x,s)>0: ‖(γμe-xs)+‖1≤βγμ},

其中:γ∈(0,1/2];β∈(0,1/2].易驗證: 如果(x,s)∈N(β,γ),則xs≥(1-β)γμe.

設(x,s)∈N(β,γ),τ∈(0,1].首先計算預估方向(Δxa,Δsa):

AΔxa-Δsa=(1-σ)τ(q0-q),sΔxa+xΔsa=(γμe-xs)-+n(γμe-xs)+,

(3)

其中σ∈(0,γ).其次計算校正方向:

AΔxc-Δsc=0,sΔxc+xΔsc=-ΔxaΔsa;

(4)

新迭代點:

其中α∈(0,1]是使得xTs和‖Ax-s+q‖能充分下降并使(x,s)∈N(β,γ)成立的最大步長.設

算法1輸入:ε>0,γ∈(0,1/2],σ∈(0,γ),β∈(0,1/4].初始點(x0,s0)∈N(β,γ),τ0=1.令k=0.當(xk)Tsk>ε時執(zhí)行如下步驟:

1) 解式(3)得到預估方向(Δxa,Δsa);

2) 解式(4)得到校正方向(Δxc,Δsc);

注1記殘余為rk=Axk-sk+q,利用式(6),(7)可得

顯然,(xk)Tsk→0?rk→0.

2 收斂性分析

引理1[2]如果LCP是P*(κ),則對任意的(x,s)>0和a,b∈n,系統(tǒng):Au-v=b,su+xv=a有唯一解(u,v),并且滿足:其中:

引理2設(x*,s*)∈F*,并且序列(x,s,τ)由算法1產(chǎn)生,則有估計:

τ((x0-x*)Ts+(s0-s*)Tx)≤(2+8κ)nμ.

證明: 顯然A(τ(x0-x*)+x*-x)-(τ(s0-s*)+s*-s)=0.再注意到xTs≥τ(x0)Ts0,利用文獻[2]中引理3.5的證明方法易證結論成立.

證明：因為Ax*-s*=q,Ax0-s0=q0,則A(x0-x*)-(s0-s*)=q0-q.因此,對于z=(x,s),有

又由于

其中最后一個不等式由引理2和(x,s)∈N(β,γ)可得.同理,有

因此結論成立.

引理4設(x,s)∈N(β,γ),(Δxa,Δsa)是式(3)的解,則存在常數(shù)C>1,使得:

證明：易知〈(γμe-xs)-,(γμe-xs)+〉=0.將引理1用于式(3),可得

記I={i:xisi≥γμ},可得

利用‖(Δxa,Δsa)‖z的定義,可知結論成立.證畢.

在上面證明中使用同一常數(shù)C≥1,后面的論述中將采用同一處理方法.

引理5設(x,s)∈N(β,γ),(Δxc,Δsc)是式(4)的解,則:

證明：由式(4)知,D-1Δxc+DΔsc=(xs)-1/2(-ΔxaΔsa).由引理1和引理4的證明,可得

推論1設(x,s)∈N(β,γ),則存在常數(shù)C≥1,使得:

1) |(Δxa)TΔsc|/n≤‖D-1Δxa‖‖DΔsc‖≤C(1+κ)5n3μ;

2) |(Δxc)TΔsa|/n≤‖D-1Δxc‖‖DΔsa‖≤C(1+κ)5n3μ;

3) |(Δxc)TΔsc|/n≤‖D-1Δxc‖‖DΔsc‖≤C(1+κ)7n4μ.

證明：僅證1),其他類似可證.首先,有

|(Δxa)TΔsc|/n=|eT(ΔxaΔsc)|/n≤‖e‖‖ΔxaΔsc‖/n≤‖ΔxaΔsc‖.

其次,由引理4和引理5,可得

‖ΔxaΔsc‖≤‖D-1Δxa‖‖DΔsc‖≤C(1+κ)5n3μ.

Δx(α)Δs(α)=α3(ΔxaΔsc+ΔsaΔxc)+α4ΔxcΔsc;μ(α)=x(α)Ts(α)/n.

引理6設(x,s)∈N(β,γ).如果α∈(0,α0],則‖Δx(α)Δs(α)‖1≤αnβγμ(α).

證明：由式(3),(4)可得x(α)s(α)=χ(α)+Δx(α)Δs(α).根據(jù)推論1,可得

其中第二個不等式用了H?lder不等式.注意到∑(γμe-xs)+≥0,于是

因此,由推論1可得

由式(8),(10),可得

其中h(α)=-βγ+αβγ+2α2C(1+κ)5n2+α3C(1+κ)7n3.對任意的α∈(0,α0],有

所以結論成立.由式(9),可得∑χ(α)≤nμ+α((γ-1)nμ+(n-1)βγμ)≤nμ+α(γ+βγ-1)nμ.

引理7設(x,s)∈N(β,γ),有μ(α)≤μ,α∈(0,α0].

證明：由式(10)中第一個等式及推論1知,對于α∈(0,α0],有

因為γ+βγ+1/5+1/100∈(0,1),所以對任意的α∈(0,α0],有μ(α)≤μ.

引理8設(x,s)∈N(β,γ).若μ(α)≤μ,則‖(γμ(α)e-χ(α))+‖1≤(1-nα)βγμ(α)e.

易知,χ(α)≥0.因此,

定理1設(x,s)∈N(β,γ),則由式(5)所定義的αc滿足αc≥α0.

證明：由引理6和引理8知,對于α∈(0,α0],有

所以有x(α)s(α)≥(1-β)γμ(α)e>0.又因為x>0,s>0,由連續(xù)性可知α∈(0,α0],從而有x(α)>0,s(α)>0.證畢.

不失一般性,下面設σ=λ/2.

定理2設(x,s)∈N(β,γ),則有αf≥α0.

證明：由式(9),(10)及推論1,有

x(α)Ts(α)-(1-(1-σ)α)xTs≥αγnμ/2-2α3C(1+κ)5n3μ-α4C(1+κ)7n4μ∶=αnμf(α),

其中f(α)=γ/2-2α2C(1+κ)5n2-α3C(1+κ)7n3.對于α∈(0,α0],有f(α)≥f(α0)≥γ(1/2-2/C-1/C2)>0.結論成立.

定理3設(x,s)∈N(β,γ),則αg≥α0.

證明：令δ=γ+βγ+1/5+1/100.由引理7可證結論成立.

定理4算法1至多需要O((1+κ)5/2nlogε-1)次迭代終止.

3 數(shù)值結果

為檢驗算法的有效性,通過編寫Matlab程序做如下數(shù)值實驗.取γ=0.005,β=0.001,σ=0.000 1,ε≤1×10-7.

表1 例1的數(shù)值結果Table 1 Results of example 1

例2[10]取A=MTM+N-NT+D,其中:M=5-10rand(n);N=5rand(n);D是對角矩陣,且Di在(0,0.3)上平均分布,q的分量在(-500,500)上平均分布.利用算法1所得數(shù)值結果列于表2.

表2 例2的數(shù)值結果Table 2 Results of example 2

綜上,本文分析了求解P*(κ)-LCP的二階預估-校正內(nèi)點算法,該算法具有目前不可行算法最好的多項式復雜度O((1+κ)5/2nL).數(shù)值實驗驗證了算法的有效性.

[1] Karmarkar N K.A New Polynomial-Time Algorithm for Linear Programming [J].Combinatorica,1984,4(4): 373-395.

[2] Potra F A,Stoer J.On a Class of Superlinearly Convergent Polynomial Time Interior Point Methods for Sufficient LCP [J].SIAM J Optim,2009,20(3): 1333-1363.

[3] Mehrotra S.On the Implementation of a Primal-Dual Interior Point Method [J].SIAM J Optim,1992,2(4): 575-601.

[6] Bai Y Q,Ghami M E,Roos C.A Comparative Study of Kernel Functions for Primal-Dual Interior-Point Algorithms in Linear Optimization [J].SIAM J Optim,2004,15(1): 101-128.

[9] LIU Hong-wei,LIU Xin-ze,LIU Chang-he.Mehrotra-Type Predictor-Corrector Algorithms for Sufficient Linear Complementarity Problem [J].Appl Numer Math,2012,62(12): 1685-1700.

[10] Kanzow C.Global Convergence Properties of Some Iterative Methods for Linear Complementarity Problems [J].SIAM J Optim,1996,6(2): 326-341.