一种修正三项Hestenes-Stiefel共轭梯度投影算法及其应用

李丹丹, 李远飞, 王松华

(1. 广州华商学院 应用数学系, 广州 511300; 2. 百色学院 数学与统计学院, 广西 百色 533000)

0 引 言

非线性方程组广泛应用于通信、 化学、 电力系统、 压缩感知等领域[1-3], 因此建立高效求解非线性方程组的数值方法具有重要意义. 本文主要考虑如下非线性方程组：

H(z)=0,z∈n,

(1)

其中H:n→n为连续单调函数, 即对于任意的z1,z2∈n, 有

(H(z1)-H(z2))T(z1-z2)≥0.

(2)

目前, 求解非线性方程组的方法主要有牛顿法、 拟牛顿法、 共轭梯度法等[4-5], 其中牛顿法和拟牛顿法由于局部收敛较快等优势, 广泛应用于求解中小规模的非线性方程组中, 但在运算过程中每次迭代都需要计算一个与Jacobi矩阵相关的线性方程组, 运算复杂性较高, 限制了其在求解大规模非线性方程组问题中的应用.共轭梯度法由于在求解过程中无需计算方程组的Jacobi矩阵、 编程简单且存储空间要求较小等特点而成为求解大规模优化问题的主要方法之一, 目前已成功地应用于求解大规模非线性方程问题和压缩感知问题中[6-8].

不失一般性, 共轭梯度法的迭代式为zk+1=zk+tkdk, 其中:zk为当前迭代点;zk+1为下一个迭代点;tk为某种线搜索产生的步长;dk为某种方法产生的搜索方向, 其迭代形式为

式中βk为共轭参数,H(zk)简记为Hk.显然, 共轭梯度法的优劣受给定线搜索方法和搜索方向的影响, 而搜索方向的好坏取决于共轭参数, 其决定了算法的收敛速度.

针对单调非线性优化问题, 文献[9]采用一种高效的线搜索技术, 令步长tk=srmk,mk为满足下列不等式的最小非负整数：

-H(zk+tkdk)Tdk≥βtk‖dk‖2,

(3)

其中s>0, 00.

经典Hestenes-Stiefel(HS)方法因自动满足共轭条件、 数值结果较好等优势而备受关注, 但研究表明, 该方法收敛性质一般. 为改善HS算法的收敛性质, 文献[10]在经典HS方法的基础上, 通过添加一个扰动项, 给出了一种三项搜索方向, 其表达式为

(4)

Rk={z∈n|H(uk)T(z-uk)=0},

将最优解z*与当前迭代点zk严格分离, 下一个新的迭代点zk+1可由下式计算得到：

(5)

基于上述分析, 本文在式(4)的基础上, 通过构建一个具有良好性质的新型搜索方向, 利用式(3)的线搜索技术和式(5)的超平面投影方法, 提出一种无导数修正三项HS共轭梯度投影方法, 分析新算法的全局收敛性, 给出数值实验结果, 并将其应用于信号恢复问题中.

1 算法设计

(6)

因此, 式(4)可改写为

(7)

(8)

其中δ>0为给定的常数.式(8)所构造的搜索方向不仅保留了文献[10]的充分下降性, 而且克服了文献[10]中定义的不适定性等问题.

通过上述讨论, 可建立一种修正三项HS共轭梯度投影算法, 用于求解非线性单调方程组(1), 算法描述如下.

算法1MHS(Modified Hestenes-Stiefel)算法.

步骤1) 给定初始点z0∈n(k=0), 运行参数s>0, 00,δ>0；

步骤2) 计算函数值H(zk), 如果‖H(zk)‖≤ε, 则算法停止, 否则转步骤3)；

步骤3) 利用式(6)和式(8)计算搜索方向dk；

步骤4) 利用式(3)的线搜索方法计算步长tk, 令uk=zk+tkdk；

步骤5) 计算函数值H(uk), 如果‖H(uk)‖≤ε, 则算法停止, 否则由式(5)计算出新的迭代点zk+1, 令k∶=k+1, 转步骤2).

2 全局收敛性分析

下面分析MHS算法的全局收敛性质.

假设:

(H1) 非线性单调方程组(1)的解集非空；

(H2) 函数H(z)在n上是Lipschitz连续的, 即存在常数a>0, 使得

‖H(z1)-H(z2)‖≤a‖z1-z2‖, ∀z1,z2∈n.

(9)

引理1若搜索方向dk由式(6)和式(8)产生, 则搜索方向dk满足充分下降性, 即对于任意的k, 均有

(10)

于是式(10)得证.

利用Cauchy-Schwatz不等式, 由式(10)易推出

‖Hk‖≤‖dk‖.

(11)

下面证明MHS算法的适定性并给出步长tk的下界.

引理2如果假设(H1),(H2)成立, 则存在一个步长tk使得式(3)成立.

证明： 类似于文献[13]中引理3.1可证, 故略.

引理3如果假设(H1),(H2)成立, 且MHS算法产生序列{zk}和{uk}, 则有

联合式(10)和式(9), 得

证明： 先证序列{zk}和{uk}有界.假设方程组(1)的最优解为z*, 即H(z*)=0, 由式(2)得(H(uk)-H(z*))T(uk-z*)≥0, 可化简为H(uk)T(zk-z*+uk-zk)≥0, 于是有

H(uk)T(zk-z*)≥H(uk)T(zk-uk).

(12)

此外, 由式(3)的线搜索方法及uk的定义得

(13)

观察式(14)可知, {‖zk-z*‖}是关于k的单调不增序列且收敛, 因此序列{zk}有界.由式(14)还可得：

‖zk+1-z*‖2≤‖zk-z*‖2,

(15)

利用递推关系可得

‖zk-z*‖2≤‖zk-1-z*‖2≤…≤‖z0-z*‖2,

结合假设(H1),(H2)可得

‖H(zk)‖=‖H(zk)-H(z*)‖≤a‖z0-z*‖.

令c=a‖z0-z*‖, 则序列{H(zk)}有界, 即对于任意的k, 均有

‖H(zk)‖≤c.

(16)

由式(13)和Cauchy-Schwatz不等式得

再联合式(16)和序列{zk}的有界性, 可知序列{uk}也有界.

因为序列{uk}有界, 所以{‖uk-z*‖}也有界, 故存在一个正常数e, 使得‖uk-z*‖≤e成立, 再结合假设(H1),(H2)得

‖H(uk)‖=‖H(uk)-H(z*)‖≤a‖uk-z*‖≤ae,

联合式(14)得

从而有

下面证明MHS算法的全局收敛性质.

证明： 用反证法.假设结论不成立, 则存在κ>0, 使得对任意的k, 均有

‖H(zk)‖≥κ.

(17)

结合式(11)得

κ≤‖Hk‖≤‖dk‖.

(18)

另一方面, 由引理4可知{zk}有界, 故存在正常数e1, 使得‖zk‖≤e1, 再结合wk-1的定义及假设(H1),(H2)得

‖wk-1‖=‖Hk-Hk-1‖≤a‖zk-zk-1‖≤a(‖zk‖+‖zk-1‖)≤2ae1.

(19)

(20)

3 数值实验

3.1 非线性单调方程组

为考察MHS算法的性能, 下面将其与MCG(Modified three-term CG method)算法[13]和MPRP(Modified Polak-Ribière-Polyak)算法[11]进行数值效果对比. 算法参数设置为s=1,r=0.5,β=0.01,δ=1.61. 测试计算机环境为Windows10(64 bite), RAM 8 GB, CPU 3.60 GHz, 用MATLAB(2014a)软件进行编程. 程序终止准则为‖Hk‖≤10-5, 迭代次数上限为3 000, 维数为4 500,12 000,24 000,30 000,45 000.测试问题选自文献[14], 函数名称和初始点列于表1.实验结果列于表2, 其中t为程序运行时间, “—”表示迭代次数已超过3 000.

表1 测试问题

表2 实验结果

续表2

由表2可见： MHS算法在规定的迭代次数内均能成功求解这10个函数, 而MCG算法和MPRP算法不能很好地求解第5个函数； MHS算法比MCG算法和MPRP算法在求解同类问题时迭代次数更少. 为更直观展示3种不同算法的优劣性, 采用文献[15]的性能曲线评价上述3种算法的性能. 图1～图3分别为不同算法的迭代次数、 函数计算次数和CPU运行时间性能比较结果, 其中曲线越靠上, 表示算法越稳定、 有效. 图1～图3中横坐标t表示给定的数值比率， 纵坐标计算公式如下：

其中P表示测试集, |P|表示测试集的问题个数,V表示算法集合,rp,s表示某种指标(如运行时间、 迭代次数和函数计算次数等).

图1 不同算法的迭代次数性能比较Fig.1 Performance comparison of number of iterations by different algorithms

图2 不同算法的函数计算次数性能比较Fig.2 Performance comparison of calculation times of functions by different algorithms

图3 不同算法的运行时间性能比较Fig.3 Performance comparison of running time by different algorithms

由图1～图3可见, MHS算法总体上比MCG算法和MPRP算法性能更优, 以更小的迭代次数和更短的运行时间求解了非线性单调方程组问题.

3.2 信号恢复问题

下面用MHS算法求解压缩感知的信号恢复问题[16], 并与MPRP算法进行比较, 以验证本文算法的性能.

表3为两种算法在求解信号恢复问题时的运行时间(t/s)、 迭代次数及MSE结果. 图4为原始信号(A)、 观测信号(B)、 用MHS算法(C)和MPRP算法(D)恢复的信号对比结果. 由表3和图4可见, MHS算法和MPRP算法均能有效地从观察信号中恢复出原始信号, 但MHS算法的求解效率更高.

表3 两种算法对稀疏信号的恢复结果

(A) 原始信号(n=4 096, 非零个数为128); (B) 观察信号; (C) MHS算法(MSE=1.317 11×10-5, 迭代次数为619, 运行时间为68.984 4 s); (D) MPRP算法(MSE=1.190 12×10-5, 迭代次数为654, 运行时间为78.937 5 s).图4 两种算法对稀疏信号的恢复对比Fig.4 Comparison of sparse signal recovery by two algorithms

实验结果表明, 本文提出的新算法不仅具有很好的理论性质, 且数值效果较好, 无论是在求解大规模非线性方程组问题还是在稀疏信号恢复问题上均有效.