賀曉瑞 湯京永

(信陽(yáng)師范大學(xué)數(shù)學(xué)與統(tǒng)計(jì)學(xué)院,信陽(yáng),464000)

1 引言

2012 年,國(guó)際著名優(yōu)化專(zhuān)家Potra 教授在文獻(xiàn)[1]中提出了加權(quán)線性互補(bǔ)問(wèn)題的數(shù)學(xué)模型:

其中P∈R(n+m)×n,Q∈R(n+m)×n,R∈R(n+m)×m為給定矩陣,a∈Rn+m為給定向量,w≥0 為權(quán)重向量,xs表示向量x和s對(duì)應(yīng)分量相乘組成的向量,即xs:=(x1s1,···,xnsn)T.如果加權(quán)線性互補(bǔ)問(wèn)題(1.1)滿足:對(duì)任意的(?x,?s,?y)∈Rn×Rn×Rm,當(dāng)P?x+Q?s+R?y=0 時(shí),總有〈?x,?s〉≥0,則稱(chēng)其是單調(diào)的.

經(jīng)濟(jì)、管理等領(lǐng)域中的許多均衡問(wèn)題都可以轉(zhuǎn)化成加權(quán)線性互補(bǔ)問(wèn)題進(jìn)行求解.雖然有些均衡問(wèn)題可以轉(zhuǎn)化成互補(bǔ)問(wèn)題,比如著名的Fisher 均衡問(wèn)題,但將其轉(zhuǎn)化成加權(quán)線性互補(bǔ)問(wèn)題可以更有效地進(jìn)行求解[1].

許多學(xué)者對(duì)加權(quán)線性互補(bǔ)問(wèn)題展開(kāi)了深入研究.Potra 在文獻(xiàn)[1]中給出了兩個(gè)求解單調(diào)加權(quán)線性互補(bǔ)問(wèn)題的內(nèi)點(diǎn)算法,分析了算法的多項(xiàng)式迭代復(fù)雜性?在文獻(xiàn)[2]中給出了一個(gè)求解充分加權(quán)線性互補(bǔ)問(wèn)題的預(yù)估校正內(nèi)點(diǎn)算法,證明了算法的迭代復(fù)雜性與求解充分線性互補(bǔ)問(wèn)題內(nèi)點(diǎn)算法的最好結(jié)果相同.Zhang[3]和Tang[4]分別研究了求解單調(diào)加權(quán)線性互補(bǔ)問(wèn)題的光滑牛頓法,分析了這些算法的全局和局部收斂性質(zhì).Asadi 等[5]給出了第一個(gè)求解單調(diào)加權(quán)線性互補(bǔ)問(wèn)題全牛頓步內(nèi)點(diǎn)算法,證明了算法具有目前最好的迭代復(fù)雜性.最近,Tang 和Zhou[6]給出了一個(gè)求解非單調(diào)加權(quán)線性互補(bǔ)問(wèn)題的阻尼高斯牛頓法,在局部誤差界條件下證明了算法具有局部二次收斂性質(zhì).

本文考慮一種特殊的加權(quán)線性互補(bǔ)問(wèn)題:

其中M∈Rn×n為給定矩陣,q∈Rn為給定向量.因?yàn)閤≥0,s≥0,xs=0 當(dāng)且僅當(dāng)x≥0,s≥0,xT s=0,所以如果權(quán)重向量w為零向量,那么互補(bǔ)問(wèn)題(1.2)即為如下標(biāo)準(zhǔn)的線性互補(bǔ)問(wèn)題:

光滑牛頓法是求解加權(quán)線性互補(bǔ)問(wèn)題(1.1)和線性互補(bǔ)問(wèn)題(1.3)的一種重要方法([3,4,7,8]),其基本思想是利用光滑函數(shù)將優(yōu)化問(wèn)題轉(zhuǎn)化成一個(gè)光滑方程組,然后利用牛頓法求解該方程組.

受文獻(xiàn)[3,4,7,8]的啟發(fā),本文利用一個(gè)帶有權(quán)重的光滑函數(shù)將加權(quán)線性互補(bǔ)問(wèn)題(1.2)轉(zhuǎn)化成一個(gè)光滑方程組,然后用一個(gè)預(yù)估校正光滑牛頓法去求解.在適當(dāng)條件下,我們將證明給出的算法具有全局和局部二次收斂性質(zhì).特別地,在解集非空的條件下,我們將證明價(jià)值函數(shù)點(diǎn)列收斂到零.最后,我們將用數(shù)值試驗(yàn)驗(yàn)證我們的算法的有效性.

2 問(wèn)題(1.2)的等價(jià)轉(zhuǎn)化

考慮如下光滑函數(shù):

其中c≥0 是一個(gè)固定常數(shù).

接下來(lái)，我們利用光滑函數(shù)(2.1)將特殊加權(quán)互補(bǔ)問(wèn)題(1.2)等價(jià)轉(zhuǎn)化成一個(gè)非線性方程組H(z)=0.

引理2.1?c在任意點(diǎn)(ε,a,b)∈R++×R×R 處連續(xù)可微且

引理2.2 對(duì)任意的(ε,a,b)∈R+×R×R,有

證明如果?c(ε,a,b)=0,則由(2.1)可得

上式兩邊同時(shí)平方,可得ab=c+ε,從而

如果a≥0,b≥0,ab=c+ε,則由(2.1)可得?c(ε,a,b)=a+b-=0.證畢.

令z:=(ε,x,s),其中ε為光滑函數(shù)中的光滑參數(shù),x,s為加權(quán)線性互補(bǔ)問(wèn)題的變量.利用(2.1)給出的函數(shù)?,我們定義

其中

由引理2.2 可知H(z)=0 當(dāng)且僅當(dāng)ε=0 并且(x,s)是特殊加權(quán)線性互補(bǔ)問(wèn)題(1.2)的解.

引理2.3 以下結(jié)論成立:

(i)H(z)在任意點(diǎn)z∈R++×Rn×Rn處連續(xù)可微,其雅可比矩陣為

這里I表示n×n單位矩陣,而

其中vec{αi}表示向量α=(α1,···,αn)T,diag{αi}表示第i個(gè)對(duì)角元素為αi的對(duì)角矩陣.

(ii)如果M半正定,那么H′(z)在任意的點(diǎn)z∈R++×Rn×Rn處非奇異.

證明由引理2.1 可知結(jié)論(i)成立.類(lèi)似于文獻(xiàn)[3]中的引理1,我們可證明結(jié)論(ii).證畢.

3 算法的構(gòu)造

下面我們給出求解問(wèn)題(1.2)的算法.

算法3.1選擇參數(shù)σ,δ∈(0,1)和初始點(diǎn)z0:=(ε0,x0,s0)∈R++×Rn×Rn?選取γ∈(0,1)使得γ <ε0且γ+σ <1?選取θ >0 和序列 {θk}使得≤θ <∞.令k:=0.

步驟1 如果∥H(zk)∥=0,則停止迭代.

步驟2(預(yù)估步) 如果∥H(zk)∥>1,則令:=zk并轉(zhuǎn)步驟3.否則,通過(guò)求解線性系統(tǒng)

令αk:=,其中l(wèi)k是滿足下式的最小非負(fù)整數(shù)l,

步驟4 令k:=k+1.轉(zhuǎn)步驟1.

注3.1 算法3.1 的框架由Huang,Han 和Chen[9]提出,用于求解線性互補(bǔ)問(wèn)題(1.3).與文獻(xiàn)[9]中算法的不同之處是,在算法3.1 的校正步中,我們采用非單調(diào)線搜索技術(shù)來(lái)生成步長(zhǎng),這使得我們的算法具有更好的計(jì)算效果.

定理3.1 如果M半正定,那么算法3.1 可以產(chǎn)生一個(gè)無(wú)窮點(diǎn)列 {zk=(εk,xk,sk)},并且對(duì)所有的k≥0 有εk >0.

因此,我們得到如下結(jié)論:如果zk∈R++×Rn×Rn,那么算法3.1 可以生成zk+1并且zk+1∈R++×Rn×Rn.因?yàn)閦0∈R++×Rn×Rn,故由數(shù)學(xué)歸納法可知定理成立.證畢.

4 收斂性分析

引理4.1 設(shè){zk=(εk,xk,sk)}是由算法3.1 生成的迭代序列,則對(duì)所有的k≥0 有：(i)∥H(zk+1)∥≤(1+θk)∥H(zk)∥;(ii)εk≥βk.

證明由(3.5)可知,對(duì)所有的k≥0 有∥H(zk+1)∥≤.如果預(yù)估步成立,那么∥H(zk)∥≤1 且

否則,

結(jié)論(i)得證.

接下來(lái),我們假設(shè)εk≥βk對(duì)某個(gè)k成立,例如,它對(duì)k=0 成立.如果預(yù)估步成立,那么由(3.1)和(3.2)可得

這里最后一個(gè)不等式成立是因?yàn)樾蛄衶βk}單調(diào)遞減.由數(shù)學(xué)歸納法可知結(jié)論(ii)成立.證畢.

定理4.1 假設(shè)M是半正定的并且{zk=(εk,xk,sk)}是由算法3.1 生成的迭代序列,那么{zk}的任意聚點(diǎn)z?:=(ε?,x?,s?)是H(z)=0 的一個(gè)解.

證明對(duì)所有的k≥0,因?yàn)椤蜨(zk+1)∥≤(1+θk)∥H(zk)∥且<∞,所以由文獻(xiàn)[10]中的引理2.2 可得,存在一個(gè)常數(shù)H?≥0 使得

由(4.1)和(4.2),對(duì)所有的k≥0 有

現(xiàn)在我們假設(shè)H(z?)0,然后推出一個(gè)矛盾.顯然,∥H(z?)∥=H?>0.由算法3.1 的步驟3 可得

由此可得

因?yàn)閧βk}是單調(diào)遞減的,所以存在一個(gè)常數(shù)β?≥0 使得=β?.因?yàn)?H?>0,由(3.4)可得β?>0.結(jié)合引理4.1(ii)可知

因此,H(z)在z?處連續(xù)可微.在(4.9)中令k(∈K)→∞,可得

另一方面,由(3.3)可知

這里不等式成立是因?yàn)棣?≤∥H(z?)∥和β?≤γ∥H(z?)∥.由(4.10) 和(4.11) 可知(1-γ-σ)∥H(z?)∥≤0,再結(jié)合γ+σ <1 可得H(z?)=0,從而推出矛盾.證畢.

下面我們證明價(jià)值函數(shù)點(diǎn)列 {∥H(zk)∥} 收斂到零.

引理4.2 設(shè)Φw由(2.3)定義.則對(duì)任意的(ε,x,s,t)∈R++×Rn×Rn×Rn,有

這里E:=(1,···,1)T.

證明由(2.1)可知對(duì)任意的(ε,a,b,?)∈R++×R3,有

從而由引理2.2 可得

利用上述結(jié)論,由(2.3)我們可以得到

證畢.

引理4.3 假設(shè)M半正定并且 {zk=(εk,xk,sk)} 是由算法3.1 生成的迭代序列.如果特殊加權(quán)線性互補(bǔ)問(wèn)題(1.2)的解集非空,那么有

證明由(4.7)可得

因?yàn)閷?duì)所有的k≥0,有∥H(zk+1)∥≤(1+θk)∥H(zk)∥,故由文獻(xiàn)[10] 中的引理2.1 可知∥H(zk)∥≤eθ∥H(z0)∥,進(jìn)而可得

這說(shuō)明序列{εk}是有界的.因?yàn)閱?wèn)題(1.2)的解集非空,所以存在(x?,s?)∈Rn×Rn使得

對(duì)任意的k=0,1,···,令

由(2.2),(4.13)和(4.14)可知,對(duì)所有的k=0,1,···,有

這說(shuō)明 {(uk,vk)} 是有界的.現(xiàn)在,我們構(gòu)造另一個(gè)序列:

由(4.17)可知Φw(εk,xk,sk)=2(εkvk-εkxk),再結(jié)合引理4.2 可得

由(4.15)可知x?s?=w,故〈x?,s?〉=.因此,

由(4.15),(4.18),(4.19),(4.23)和M半正定,可得

這里

由(4.13)和(4.24),對(duì)所有的k≥0,有

因?yàn)樾蛄?{(εk,uk,vk)} 是有界的,所以如果∥xk∥→∞,那么是有界的.此時(shí),當(dāng)k→∞時(shí)(4.26) 的左邊趨向于+∞,這與(4.26) 矛盾.因此,序列 {xk} 是有界的.再結(jié)合 {(εk,uk)} 的有界性和sk=Mxk+q+εkuk可得 {sk} 是有界的.因此,序列 {zk=(εk,xk,sk)} 有界,從而可知 {zk} 至少有一個(gè)聚點(diǎn)z?.由H的連續(xù)性和(4.12)可知H(z?)=H?.再結(jié)合定理4.1 可得H?=0,即=0,進(jìn)而由(3.4)可知β?=0,從而推出矛盾.證畢.

因?yàn)樗惴?.1 的局部二次收斂速度取決于預(yù)估步,所以由文獻(xiàn)[9]中的定理5.1 可得如下結(jié)論.

引理4.4 假設(shè)M半正定并且 {zk} 是由算法3.1 產(chǎn)生的迭代序列.令z?為 {zk} 的任意聚點(diǎn).如果所有的V∈?H(z?)都是非奇異的,則 {zk} 收斂于z?并且∥zk+1-z?∥=O(∥zk-z?∥2).

5 數(shù)值結(jié)果

在本節(jié)中,我們通過(guò)數(shù)值試驗(yàn)來(lái)檢驗(yàn)算法的有效性.算法3.1 中的參數(shù)設(shè)置為δ=0.5,ε0=10-4,γ=10-5,σ=10-5,θk=0.9k.我們使用 ∥H(zk)∥≤10-6作為停止準(zhǔn)則.

5.1 求解特殊加權(quán)線性互補(bǔ)問(wèn)題(1.2)

本小節(jié)應(yīng)用算法3.1 求解特殊加權(quán)線性互補(bǔ)問(wèn)題(1.2).令M=UUT/∥UUT∥,其中U=rand(n,n),q=rand(n,1).此外,選擇=rand(n,1),=rand(n,1),并令w:=.我們分別使用以下兩個(gè)初始點(diǎn)求解這些測(cè)試問(wèn)題：

SP1x0=s0=(1,···,1)T?

SP2x0=(1,0,···,0)T,s0=Mx0+q.

首先,我們生成一個(gè)n=1000 的測(cè)試問(wèn)題,然后應(yīng)用算法3.1 去求解此問(wèn)題.表1 給出了迭代過(guò)程中∥H(zk)∥的值.由表1,我們可以很容易看出算法3.1 具有局部快速收斂性質(zhì).

表1 第k 次迭代時(shí)∥H(zk)∥的值

接下來(lái),對(duì)于每個(gè)測(cè)試問(wèn)題,我們生成10 個(gè)算例,然后應(yīng)用算法3.1 去求解這些算例.數(shù)值試驗(yàn)的結(jié)果列于表2,其中SP 表示初始點(diǎn),AIT 和ACPU分別表示用算法3.1 求解問(wèn)題所需的迭代次數(shù)和CPU 時(shí)間的平均值,AHK 表示算法終止時(shí) ∥H(zk)∥的平均值.從表2 可以看出,算法3.1 對(duì)于求解特殊加權(quán)線性互補(bǔ)問(wèn)題(1.2)是非常有效的,它只需很少的迭代次數(shù)和CPU 時(shí)間就可以找到滿足終止條件的解.此外,我們還可以發(fā)現(xiàn),當(dāng)測(cè)試問(wèn)題的規(guī)模變大時(shí),算法3.1 所需的迭代次數(shù)變化不大，但所需的CPU 時(shí)間增加.

表2 求解問(wèn)題(1.2)的數(shù)值結(jié)果

5.2 求解加權(quán)線性互補(bǔ)問(wèn)題(1.1)

本小節(jié)應(yīng)用算法3.1 求解加權(quán)線性互補(bǔ)問(wèn)題(1.1),其中

其中N是一個(gè)給定的n×n對(duì)稱(chēng)半正定矩陣,A∈Rm×n是一個(gè)m

在試驗(yàn)中,我們產(chǎn)生一個(gè)行滿秩的矩陣A∈Rm×n,然后選擇N=BBT/∥BBT∥,其中B=rand(n,n),=rand(n,1),f=rand(n,1),再定義.對(duì)于每個(gè)問(wèn)題,我們分別使用以下兩個(gè)初始點(diǎn)進(jìn)行求解:

SP1x0=s0=(1,0,···,0)T,y0=(0,···,0)T?

SP2x0=rand(n,1),s0=rand(n,1),y0=rand(m,1).

首先,我們生成一個(gè)n=1000,m=500 的測(cè)試問(wèn)題,然后應(yīng)用算法3.1 去求解.表3 給出了迭代過(guò)程中∥H(zk)∥的值,其再一次表明算法3.1 具有局部快速收斂速度.

表3 第k 次迭代時(shí)∥H(zk)∥的值

接下來(lái),我們對(duì)算法3.1 和文獻(xiàn)[9],[11]給出的算法進(jìn)行比較.

本文給出的算法3.1,記為NPCM?文獻(xiàn)[9]給出的求解線性互補(bǔ)問(wèn)題的預(yù)估校正光滑牛頓法,記為PCM?文獻(xiàn)[11]給出的求解加權(quán)互補(bǔ)問(wèn)題的非單調(diào)光滑牛頓法,記為NSNM.

對(duì)于每個(gè)n(=2m)的測(cè)試問(wèn)題,我們生成10 個(gè)算例,并分別應(yīng)用上述三種算法求解.數(shù)值結(jié)果列于表4.

表4 NPCM,PCM和NSNM的數(shù)值比較結(jié)果

由表4,我們有如下結(jié)論：

(i)算法3.1 對(duì)于求解一般的加權(quán)線性互補(bǔ)問(wèn)題(1.1)也是有效的,它只需很少的迭代次數(shù)和CPU 時(shí)間就可找到滿足精度要求的解.

(ii)算法3.1 比文獻(xiàn)[9]中的預(yù)估校正光滑牛頓法更穩(wěn)定.

(iii)與文獻(xiàn)[11]中的算法相比,算法3.1 所需的迭代次數(shù)減少而CPU 時(shí)間增加.一個(gè)合理的解釋是算法3.1 在每次迭代時(shí)需要求解兩個(gè)線性方程組并進(jìn)行兩次線搜索,而文獻(xiàn)[11]中的算法在每次迭代時(shí)只需求解一個(gè)線性方程組并進(jìn)行一次線搜索.換句話說(shuō),與文獻(xiàn)[11]中的算法相比,算法3.1 每一次迭代時(shí)∥H(z)∥下降的量更大,但所需的CPU 時(shí)間更多.