李丹丹, 李遠飛, 王松華
(1. 廣州華商學院 應用數學系, 廣州 511300; 2. 百色學院 數學與統(tǒng)計學院, 廣西 百色 533000)
非線性方程組廣泛應用于通信、 化學、 電力系統(tǒng)、 壓縮感知等領域[1-3], 因此建立高效求解非線性方程組的數值方法具有重要意義. 本文主要考慮如下非線性方程組:
H(z)=0,z∈n,
(1)
其中H:n→n為連續(xù)單調函數, 即對于任意的z1,z2∈n, 有
(H(z1)-H(z2))T(z1-z2)≥0.
(2)
目前, 求解非線性方程組的方法主要有牛頓法、 擬牛頓法、 共軛梯度法等[4-5], 其中牛頓法和擬牛頓法由于局部收斂較快等優(yōu)勢, 廣泛應用于求解中小規(guī)模的非線性方程組中, 但在運算過程中每次迭代都需要計算一個與Jacobi矩陣相關的線性方程組, 運算復雜性較高, 限制了其在求解大規(guī)模非線性方程組問題中的應用.共軛梯度法由于在求解過程中無需計算方程組的Jacobi矩陣、 編程簡單且存儲空間要求較小等特點而成為求解大規(guī)模優(yōu)化問題的主要方法之一, 目前已成功地應用于求解大規(guī)模非線性方程問題和壓縮感知問題中[6-8].
不失一般性, 共軛梯度法的迭代式為zk+1=zk+tkdk, 其中:zk為當前迭代點;zk+1為下一個迭代點;tk為某種線搜索產生的步長;dk為某種方法產生的搜索方向, 其迭代形式為
式中βk為共軛參數,H(zk)簡記為Hk.顯然, 共軛梯度法的優(yōu)劣受給定線搜索方法和搜索方向的影響, 而搜索方向的好壞取決于共軛參數, 其決定了算法的收斂速度.
針對單調非線性優(yōu)化問題, 文獻[9]采用一種高效的線搜索技術, 令步長tk=srmk,mk為滿足下列不等式的最小非負整數:
-H(zk+tkdk)Tdk≥βtk‖dk‖2,
(3)
其中s>0, 0
經典Hestenes-Stiefel(HS)方法因自動滿足共軛條件、 數值結果較好等優(yōu)勢而備受關注, 但研究表明, 該方法收斂性質一般. 為改善HS算法的收斂性質, 文獻[10]在經典HS方法的基礎上, 通過添加一個擾動項, 給出了一種三項搜索方向, 其表達式為
(4)
Rk={z∈n|H(uk)T(z-uk)=0},
將最優(yōu)解z*與當前迭代點zk嚴格分離, 下一個新的迭代點zk+1可由下式計算得到:
(5)
基于上述分析, 本文在式(4)的基礎上, 通過構建一個具有良好性質的新型搜索方向, 利用式(3)的線搜索技術和式(5)的超平面投影方法, 提出一種無導數修正三項HS共軛梯度投影方法, 分析新算法的全局收斂性, 給出數值實驗結果, 并將其應用于信號恢復問題中.
(6)
因此, 式(4)可改寫為
(7)
(8)
其中δ>0為給定的常數.式(8)所構造的搜索方向不僅保留了文獻[10]的充分下降性, 而且克服了文獻[10]中定義的不適定性等問題.
通過上述討論, 可建立一種修正三項HS共軛梯度投影算法, 用于求解非線性單調方程組(1), 算法描述如下.
算法1MHS(Modified Hestenes-Stiefel)算法.
步驟1) 給定初始點z0∈n(k=0), 運行參數s>0, 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).
下面分析MHS算法的全局收斂性質.
假設:
(H1) 非線性單調方程組(1)的解集非空;
(H2) 函數H(z)在n上是Lipschitz連續(xù)的, 即存在常數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)的最優(yōu)解為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)
為考察MHS算法的性能, 下面將其與MCG(Modified three-term CG method)算法[13]和MPRP(Modified Polak-Ribière-Polyak)算法[11]進行數值效果對比. 算法參數設置為s=1,r=0.5,β=0.01,δ=1.61. 測試計算機環(huán)境為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 實驗結果
續(xù)表2
由表2可見: MHS算法在規(guī)定的迭代次數內均能成功求解這10個函數, 而MCG算法和MPRP算法不能很好地求解第5個函數; MHS算法比MCG算法和MPRP算法在求解同類問題時迭代次數更少. 為更直觀展示3種不同算法的優(yōu)劣性, 采用文獻[15]的性能曲線評價上述3種算法的性能. 圖1~圖3分別為不同算法的迭代次數、 函數計算次數和CPU運行時間性能比較結果, 其中曲線越靠上, 表示算法越穩(wěn)定、 有效. 圖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算法性能更優(yōu), 以更小的迭代次數和更短的運行時間求解了非線性單調方程組問題.
下面用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
實驗結果表明, 本文提出的新算法不僅具有很好的理論性質, 且數值效果較好, 無論是在求解大規(guī)模非線性方程組問題還是在稀疏信號恢復問題上均有效.