一个求解特殊加权线性互补问题的预估校正光滑牛顿法

贺晓瑞 汤京永

(信阳师范大学数学与统计学院,信阳,464000)

1 引言

2012 年,国际著名优化专家Potra 教授在文献[1]中提出了加权线性互补问题的数学模型:

其中P∈R(n+m)×n,Q∈R(n+m)×n,R∈R(n+m)×m为给定矩阵,a∈Rn+m为给定向量,w≥0 为权重向量,xs表示向量x和s对应分量相乘组成的向量,即xs:=(x1s1,···,xnsn)T.如果加权线性互补问题(1.1)满足:对任意的(∆x,∆s,∆y)∈Rn×Rn×Rm,当P∆x+Q∆s+R∆y=0 时,总有〈∆x,∆s〉≥0,则称其是单调的.

经济、管理等领域中的许多均衡问题都可以转化成加权线性互补问题进行求解.虽然有些均衡问题可以转化成互补问题,比如著名的Fisher 均衡问题,但将其转化成加权线性互补问题可以更有效地进行求解[1].

许多学者对加权线性互补问题展开了深入研究.Potra 在文献[1]中给出了两个求解单调加权线性互补问题的内点算法,分析了算法的多项式迭代复杂性;在文献[2]中给出了一个求解充分加权线性互补问题的预估校正内点算法,证明了算法的迭代复杂性与求解充分线性互补问题内点算法的最好结果相同.Zhang[3]和Tang[4]分别研究了求解单调加权线性互补问题的光滑牛顿法,分析了这些算法的全局和局部收敛性质.Asadi 等[5]给出了第一个求解单调加权线性互补问题全牛顿步内点算法,证明了算法具有目前最好的迭代复杂性.最近,Tang 和Zhou[6]给出了一个求解非单调加权线性互补问题的阻尼高斯牛顿法,在局部误差界条件下证明了算法具有局部二次收敛性质.

本文考虑一种特殊的加权线性互补问题:

其中M∈Rn×n为给定矩阵,q∈Rn为给定向量.因为x≥0,s≥0,xs=0 当且仅当x≥0,s≥0,xT s=0,所以如果权重向量w为零向量,那么互补问题(1.2)即为如下标准的线性互补问题:

光滑牛顿法是求解加权线性互补问题(1.1)和线性互补问题(1.3)的一种重要方法([3,4,7,8]),其基本思想是利用光滑函数将优化问题转化成一个光滑方程组,然后利用牛顿法求解该方程组.

受文献[3,4,7,8]的启发,本文利用一个带有权重的光滑函数将加权线性互补问题(1.2)转化成一个光滑方程组,然后用一个预估校正光滑牛顿法去求解.在适当条件下,我们将证明给出的算法具有全局和局部二次收敛性质.特别地,在解集非空的条件下,我们将证明价值函数点列收敛到零.最后,我们将用数值试验验证我们的算法的有效性.

2 问题(1.2)的等价转化

考虑如下光滑函数:

其中c≥0 是一个固定常数.

接下来，我们利用光滑函数(2.1)将特殊加权互补问题(1.2)等价转化成一个非线性方程组H(z)=0.

引理2.1ϕc在任意点(ε,a,b)∈R++×R×R 处连续可微且

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

证明如果ϕc(ε,a,b)=0,则由(2.1)可得

上式两边同时平方,可得ab=c+ε,从而

如果a≥0,b≥0,ab=c+ε,则由(2.1)可得ϕc(ε,a,b)=a+b-=0.证毕.

令z:=(ε,x,s),其中ε为光滑函数中的光滑参数,x,s为加权线性互补问题的变量.利用(2.1)给出的函数ϕ,我们定义

其中

由引理2.2 可知H(z)=0 当且仅当ε=0 并且(x,s)是特殊加权线性互补问题(1.2)的解.

引理2.3 以下结论成立:

(i)H(z)在任意点z∈R++×Rn×Rn处连续可微,其雅可比矩阵为

这里I表示n×n单位矩阵,而

其中vec{αi}表示向量α=(α1,···,αn)T,diag{αi}表示第i个对角元素为αi的对角矩阵.

(ii)如果M半正定,那么H′(z)在任意的点z∈R++×Rn×Rn处非奇异.

证明由引理2.1 可知结论(i)成立.类似于文献[3]中的引理1,我们可证明结论(ii).证毕.

3 算法的构造

下面我们给出求解问题(1.2)的算法.

算法3.1选择参数σ,δ∈(0,1)和初始点z0:=(ε0,x0,s0)∈R++×Rn×Rn;选取γ∈(0,1)使得γ <ε0且γ+σ <1;选取θ >0 和序列 {θk}使得≤θ <∞.令k:=0.

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

步骤2(预估步) 如果∥H(zk)∥>1,则令:=zk并转步骤3.否则,通过求解线性系统

令αk:=,其中lk是满足下式的最小非负整数l,

步骤4 令k:=k+1.转步骤1.

注3.1 算法3.1 的框架由Huang,Han 和Chen[9]提出,用于求解线性互补问题(1.3).与文献[9]中算法的不同之处是,在算法3.1 的校正步中,我们采用非单调线搜索技术来生成步长,这使得我们的算法具有更好的计算效果.

定理3.1 如果M半正定,那么算法3.1 可以产生一个无穷点列 {zk=(εk,xk,sk)},并且对所有的k≥0 有εk >0.

因此,我们得到如下结论:如果zk∈R++×Rn×Rn,那么算法3.1 可以生成zk+1并且zk+1∈R++×Rn×Rn.因为z0∈R++×Rn×Rn,故由数学归纳法可知定理成立.证毕.

4 收敛性分析

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

证明由(3.5)可知,对所有的k≥0 有∥H(zk+1)∥≤.如果预估步成立,那么∥H(zk)∥≤1 且

否则,

结论(i)得证.

接下来,我们假设εk≥βk对某个k成立,例如,它对k=0 成立.如果预估步成立,那么由(3.1)和(3.2)可得

这里最后一个不等式成立是因为序列{βk}单调递减.由数学归纳法可知结论(ii)成立.证毕.

定理4.1 假设M是半正定的并且{zk=(εk,xk,sk)}是由算法3.1 生成的迭代序列,那么{zk}的任意聚点z∗:=(ε∗,x∗,s∗)是H(z)=0 的一个解.

证明对所有的k≥0,因为∥H(zk+1)∥≤(1+θk)∥H(zk)∥且<∞,所以由文献[10]中的引理2.2 可得,存在一个常数H∗≥0 使得

由(4.1)和(4.2),对所有的k≥0 有

现在我们假设H(z∗)0,然后推出一个矛盾.显然,∥H(z∗)∥=H∗>0.由算法3.1 的步骤3 可得

由此可得

因为{βk}是单调递减的,所以存在一个常数β∗≥0 使得=β∗.因为=H∗>0,由(3.4)可得β∗>0.结合引理4.1(ii)可知

因此,H(z)在z∗处连续可微.在(4.9)中令k(∈K)→∞,可得

另一方面,由(3.3)可知

这里不等式成立是因为ε∗≤∥H(z∗)∥和β∗≤γ∥H(z∗)∥.由(4.10) 和(4.11) 可知(1-γ-σ)∥H(z∗)∥≤0,再结合γ+σ <1 可得H(z∗)=0,从而推出矛盾.证毕.

下面我们证明价值函数点列 {∥H(zk)∥} 收敛到零.

引理4.2 设Φw由(2.3)定义.则对任意的(ε,x,s,t)∈R++×Rn×Rn×Rn,有

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

证明由(2.1)可知对任意的(ε,a,b,∆)∈R++×R3,有

从而由引理2.2 可得

利用上述结论,由(2.3)我们可以得到

证毕.

引理4.3 假设M半正定并且 {zk=(εk,xk,sk)} 是由算法3.1 生成的迭代序列.如果特殊加权线性互补问题(1.2)的解集非空,那么有

证明由(4.7)可得

因为对所有的k≥0,有∥H(zk+1)∥≤(1+θk)∥H(zk)∥,故由文献[10] 中的引理2.1 可知∥H(zk)∥≤eθ∥H(z0)∥,进而可得

这说明序列{εk}是有界的.因为问题(1.2)的解集非空,所以存在(x∗,s∗)∈Rn×Rn使得

对任意的k=0,1,···,令

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

这说明 {(uk,vk)} 是有界的.现在,我们构造另一个序列:

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

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

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

这里

由(4.13)和(4.24),对所有的k≥0,有

因为序列 {(εk,uk,vk)} 是有界的,所以如果∥xk∥→∞,那么是有界的.此时,当k→∞时(4.26) 的左边趋向于+∞,这与(4.26) 矛盾.因此,序列 {xk} 是有界的.再结合 {(εk,uk)} 的有界性和sk=Mxk+q+εkuk可得 {sk} 是有界的.因此,序列 {zk=(εk,xk,sk)} 有界,从而可知 {zk} 至少有一个聚点z∗.由H的连续性和(4.12)可知H(z∗)=H∗.再结合定理4.1 可得H∗=0,即=0,进而由(3.4)可知β∗=0,从而推出矛盾.证毕.

因为算法3.1 的局部二次收敛速度取决于预估步,所以由文献[9]中的定理5.1 可得如下结论.

引理4.4 假设M半正定并且 {zk} 是由算法3.1 产生的迭代序列.令z∗为 {zk} 的任意聚点.如果所有的V∈∂H(z∗)都是非奇异的,则 {zk} 收敛于z∗并且∥zk+1-z∗∥=O(∥zk-z∗∥2).

5 数值结果

在本节中,我们通过数值试验来检验算法的有效性.算法3.1 中的参数设置为δ=0.5,ε0=10-4,γ=10-5,σ=10-5,θk=0.9k.我们使用 ∥H(zk)∥≤10-6作为停止准则.

5.1 求解特殊加权线性互补问题(1.2)

本小节应用算法3.1 求解特殊加权线性互补问题(1.2).令M=UUT/∥UUT∥,其中U=rand(n,n),q=rand(n,1).此外,选择=rand(n,1),=rand(n,1),并令w:=.我们分别使用以下两个初始点求解这些测试问题：

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

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

首先,我们生成一个n=1000 的测试问题,然后应用算法3.1 去求解此问题.表1 给出了迭代过程中∥H(zk)∥的值.由表1,我们可以很容易看出算法3.1 具有局部快速收敛性质.

表1 第k 次迭代时∥H(zk)∥的值

接下来,对于每个测试问题,我们生成10 个算例,然后应用算法3.1 去求解这些算例.数值试验的结果列于表2,其中SP 表示初始点,AIT 和ACPU分别表示用算法3.1 求解问题所需的迭代次数和CPU 时间的平均值,AHK 表示算法终止时 ∥H(zk)∥的平均值.从表2 可以看出,算法3.1 对于求解特殊加权线性互补问题(1.2)是非常有效的,它只需很少的迭代次数和CPU 时间就可以找到满足终止条件的解.此外,我们还可以发现,当测试问题的规模变大时,算法3.1 所需的迭代次数变化不大，但所需的CPU 时间增加.

表2 求解问题(1.2)的数值结果

5.2 求解加权线性互补问题(1.1)

本小节应用算法3.1 求解加权线性互补问题(1.1),其中

其中N是一个给定的n×n对称半正定矩阵,A∈Rm×n是一个m

在试验中,我们产生一个行满秩的矩阵A∈Rm×n,然后选择N=BBT/∥BBT∥,其中B=rand(n,n),=rand(n,1),f=rand(n,1),再定义.对于每个问题,我们分别使用以下两个初始点进行求解:

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

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

首先,我们生成一个n=1000,m=500 的测试问题,然后应用算法3.1 去求解.表3 给出了迭代过程中∥H(zk)∥的值,其再一次表明算法3.1 具有局部快速收敛速度.

表3 第k 次迭代时∥H(zk)∥的值

接下来,我们对算法3.1 和文献[9],[11]给出的算法进行比较.

本文给出的算法3.1,记为NPCM;文献[9]给出的求解线性互补问题的预估校正光滑牛顿法,记为PCM;文献[11]给出的求解加权互补问题的非单调光滑牛顿法,记为NSNM.

对于每个n(=2m)的测试问题,我们生成10 个算例,并分别应用上述三种算法求解.数值结果列于表4.

表4 NPCM,PCM和NSNM的数值比较结果

由表4,我们有如下结论：

(i)算法3.1 对于求解一般的加权线性互补问题(1.1)也是有效的,它只需很少的迭代次数和CPU 时间就可找到满足精度要求的解.

(ii)算法3.1 比文献[9]中的预估校正光滑牛顿法更稳定.

(iii)与文献[11]中的算法相比,算法3.1 所需的迭代次数减少而CPU 时间增加.一个合理的解释是算法3.1 在每次迭代时需要求解两个线性方程组并进行两次线搜索,而文献[11]中的算法在每次迭代时只需求解一个线性方程组并进行一次线搜索.换句话说,与文献[11]中的算法相比,算法3.1 每一次迭代时∥H(z)∥下降的量更大,但所需的CPU 时间更多.