一種求解非線性方程組的有效Levenberg-Marquardt算法

韓揚 芮紹平

摘要：通過修改LM參數，并結合非單調技術和信賴域技術給出一種求解非線性方程組的有效Levenberg-Marquardt算法，即AMLM算法。在局部誤差界條件下，證明了AMLM算法具有局部快速收斂性。數值實驗結果表明，AMLM算法穩定、有效。

關鍵詞：LM參數；非單調技術；信賴域技術；Levenberg-Marquardt算法

中圖分類號：O221.1文獻標志碼：A

非線性方程組在光學、電學、經濟、新能源、工程技術、運籌優化等方面有著廣泛的應用[1-3]。本文主要討論以下非線性方程組

f1（x1，x2， … ，xn）=0f2（x1，x2， … ，xn）=0…fm（x1，x2， … ，xn）=0

向量形式為

F（x）=0（1）

其中，x=（x1，x2， …，xn）T，F（x）：Rn→Rm連續可微，式（1）的解集記為X*，現假設X*非空。

Levenberg-Marquardt（LM）算法是求解非線性方程組的一種重要方法，在第k次迭代時，計算試探步

dk=－（JTkJk+λkI）－1JTkFk（2）

其中，Fk=F（xk），Jk=F′（xk）是F（x）在xk處的雅可比矩陣，T表示轉置，I是單位矩陣，λk≥0是LM參數，用來克服當Jk奇異或接近奇異時試探步無法求解的困難。若式（1）對應的雅可比矩陣J（x）是Lipschitz連續的且非奇異的，初始點x0接近式（1）的解x*，并且在λk選取合適的情況下，LM算法具有二次收斂性。

參數λk的選取對于LM算法是至關重要的[4-10]。文獻[4]中取λk=‖Fk‖2，在局部誤差界條件下證明了LM算法具有局部二次收斂性，在本文中，‖·‖表示2-范數；當LM算法產生的序列xk靠近解集時，λk=‖Fk‖2因太小而失去作用，當LM算法產生的序列xk遠離解集時，λk=‖Fk‖2可能很大，LM試探步小，導致LM算法全局收斂速度慢，為了減小這種影響，取λk=μk‖Fk‖[5]，其中μk通過信賴域技術逐步更新，數值實驗結果得到一定的改善。

考慮序列xk遠離解集時，LM算法的計算效率可能低，為了提高計算效率，構造

λk=μkη‖Fk‖δ1+η‖Fk‖δ（3）

其中，1≤δ≤2，η>0為正常數，μk每步由信賴域技術更新，當序列xk靠近解集時，‖Fk‖趨近于0，λk接近于μkη‖Fk‖δ，避免了λk太小而失去作用；當序列xk遠離解集，‖Fk‖可能較大，這時λk接近于μk，避免了試探步過小，進而提高計算效率。

雅可比矩陣的計算量對于LM算法的計算效率是重要的[11-14]，文獻[11]中提出改進的LM算法，在每次迭代中不僅計算一個LM步dk，還計算一個近似LM步d︿k，數值結果表明，該方法可以節省大量的雅可比矩陣計算量；為了節省更多的雅可比矩陣計算量，文獻[13]在文獻[11]的基礎上增加了一個近似LM步d～k。近些年，研究發現，使用非單調策略的算法效率更高[11-17]。綜上，本文結合非單調技術修正LM算法，并在每次迭代中計算LM步dk和近似LM步d︿k，d～k。

1 基于非單調技術和λk修正的LM算法

令φ（x）=‖F（x）‖2為式（1）的價值函數，在文獻[14]中，φ（x）第k次迭代時的實際下降量為

Aredk=‖Fk‖2－‖F（xk+dk+αkd︿k+βkd～k）‖2

預估下降量為

Predk=‖Fk‖2－‖Fk+Jkdk‖2+‖F（yk）‖2－‖F（yk）+αkJkd︿k‖2+‖F（zk）‖2－‖F（zk）+βkJkd～k‖2

比率為rk=AredkPredk。

基于非單調技術修正Aredk，并取新的實際下降量為：redk=Fl（k）2－‖F（xk+dk+αkd︿k+βkd～k）‖2，其中，Fl（k）=max0≤j≤n（k）‖Fk－j‖，k=0，1，2，…，n；n（k）=minN0，k，N0為正整數。處理后，每次迭代時，用新的比率r-k=redkPredk代替原比率rk。下面給出基于非單調技術和λk修正的LM算法，簡記為AMLM算法。

AMLM算法：

Step 1 給定x0∈Rn， ε≥0， μ0>m>0， 0

Step 2 如果‖JΤkFk‖≤ε，則停止計算；否則，轉Step 3；

Step 3 計算dk=－（JTkJk+λkI）－1JTkFk，其中λk取式（3），令yk=xk+dk，計算

d︿k=－（JTkJk+λkI）－1JTkF（yk）

令zk=xk+αkd︿k，其中αk由文獻[13]中式（22）獲得，計算d～k=－（JTkJk+λkI）－1JTkF（zk），令sk=xk+αkd︿k+βkd～k，其中βk由文獻[13]中式（28）獲得；

Step 4 計算r-k，令xk+1=xk+sk，r-k≥p0xk，??? r-k

Step 5 計算μk+1=4μk，??? r-kp2

Step 6 令k：=k+1，轉Step 2。

在AMLM算法中，要求μk不小于正常數m，防止迭代序列接近解時致使試探步過大。

2 收斂性分析

定義1 設NRn，且N∩X*≠，若存在常數c>0，使得

‖F（x）‖≥c·dist（x，X*） ，x∈N（4）

其中，dist（x，X）=infy∈X‖y－x‖，則稱F（x）在N上滿足局部誤差界條件[16]。

下文記x-k∈X，滿足‖x-k－xk‖=dist（xk，X*）。

假設1 （a）設F（x）連續可微，‖F（x）‖為式（1）在鄰域N（x*，b）上的一個局部誤差界，其中N（x*，b）=x∈Rn| ‖x－x*‖≤b， 0

（b）F（x），J（x）在N（x*，b）上是Lipcshitz連續的，存在L1>0，L2>0使得

‖J（y）－J（x）‖≤L1‖y－x‖，? x，y∈N（x*，b）（5）

‖F（y）－F（x）‖≤L2‖y－x‖，x，y∈N（x*，b）（6）

進而有

‖F（y）－F（x）－J（x）（y－x）‖≤L1‖y－x‖2， x，y∈N（x*，b）（7）

引理1 若假設1成立，xk∈N（x*，b），對所有充分大的k，有

（a）存在一個常數M>m>0，使得μk≤M；

（b）mγ≤λk≤ML2‖x-k－xk‖δ，其中γ=minη1+η，ηcδ1+η‖xk－x-k‖δ。

證明：先證（a），由文獻[13]中引理2可知，當k充分大，有

|rk－1|=Aredk－PredkPredk≤O（‖x-k－xk‖·‖dk‖2）O（‖x-k－xk‖·‖dk‖）→0

即rk→1，故

r︿k=redkPredk=F2l（k）－‖F（xk+dk+αkd︿k+βkd～k）‖2Predk≥‖Fk‖2－‖Fk+1‖2Predk=rk→1

因此，對所有充分大的k，存在一個常數M>m>0使μk≤M成立。

下面證明（b），由式（6），引理1（a）及‖F（xk）‖=0，對所有充分大的k，有

λk=μkη‖Fk‖δ1+η‖Fk‖δ≤μkη‖Fk‖δ≤MηLδ2‖xk－x-k‖δ

若‖Fk‖≤1，則‖Fk‖δ≤1，1+η‖Fk‖δ≤1+η，由式（4），得

λk=μkη‖Fk‖δ1+η‖Fk‖δ≥μkη1+η‖Fk‖δ≥mηcδ1+η‖xk－x-k‖δ

若‖Fk‖>1，則‖Fk‖δ>1，1+η‖Fk‖δ≤（1+η）‖Fk‖δ，故

λk=μkη‖Fk‖δ1+η‖Fk‖δ≥μkη‖Fk‖δ（1+η）‖Fk‖δ≥mη1+η

令γ=minη1+η，ηcδ1+η‖xk－x-k‖δ，得λk≥mγ。

綜上，對所有充分大的k，有mγ≤λk≤ML2‖x-k－xk‖δ。

引理2 若假設1成立，且k充分大，則存在正常數c1，c2，c3使得

（a） ‖dk‖≤c1dist（xk，X*）；（b）‖d︿k‖≤c2dist（xk，X*）；（c）‖d～k‖≤c3dist（xk，X*）

證明：先證（a），令φk（d）=‖Fk+Jkd‖2+λk‖d‖2，由式（2）知，AMLM算法產生的dk是φk（d）的極小點，再由式（4），式（6）及式（7）得

‖dk‖2≤1λkφk（xk－x-k）

=1λk（‖Fk+Jk（xk－x-k）‖2+λk‖xk－x-k‖2）

=1+η‖Fk‖δμkη‖Fk‖δ‖Fk+Jk（xk－x-k）‖2+‖xk－x-k‖2

≤1+ηLδ2‖xk－x-k‖δmηcδ‖xk－x-k‖δL12‖xk－x-k‖4+‖xk－x-k‖2

故存在正常數c1使‖dk‖≤c1dist（xk，X*）。

下證（b）和（c），首先研究‖（JTkJk+λkI）－1JTk‖，設J（x-）的奇異值分解為

J（x-）=U-Σ-V-=U-1，U-2Σ-1

0V-1V-2=U-1Σ-1V-1

其中，Σ-1=diag（σ-1，σ-2，…，σ-r），σ-1≥σ-2≥…≥σ-r>0，rank（J（x-））=r。相應地，對J（x）

J（x）=UΣVT=U1，U2，U3Σ1Σ2Σ3V-1V-2V-3=U1Σ1VT1+U2Σ2VT2（8）

其中，Σ1=diag（σ1，σ2，…，σr），Σ2=diag（σr+1，σr+2，…，σr+q），σ1≥σ2≥…≥σr>0，σr+1≥σr+2≥…≥σr+q>0。

因此

‖（JTkJk+λkI）－1JTk‖≤ ‖（Σ1+λkI）－1Σ1‖+‖λk－1Σ2‖（9）

根據矩陣攝動定理[18]及J（x）是Lipschitz連續，有

‖diag（Σ1－Σ-1，Σ2，0）‖≤ ‖Jk－J-k‖≤L1‖xk－x-k‖

則

‖Σ1－Σ-1‖≤L1‖xk－x-k‖，‖Σ2‖≤ ‖Jk－J-k‖≤L1‖xk－x-k‖（10）

故

‖λ－1kΣ2‖=‖Σ2‖‖λk‖≤L1mγ‖x-k－xk‖（11）

由式（10），有

‖Σ－11‖≤1σ-r－L1‖xk－x-k‖≤2σ-r（12）

對任意的σi（i=1，2，…，r），σiσ2i+λk≤σi2σiλk=12λk，那么‖（Σ21+λkI）－1Σ1‖≤12mγ，結合式（9）和（11），有

‖（JTkJk+λkI）－1JTk‖≤‖（Σ21+λkI）－1Σ1‖+‖λ－1kΣ2‖≤12mγ+L1mγ‖x-k－xk‖（13）

由d︿k的定義及式（7），得

‖d︿k‖=‖－（JTkJk+λkI）－1JTkF（yk）‖≤2‖dk‖+L1‖dk‖2‖（JTkJk+λkI）－1JTk‖

再由（a）和式（13），存在正常數c2使‖d︿k‖≤c2dist（xk，X*），同理，存在正常數c3使得

‖d～k‖=‖－（JTkJk+λkI）－1JTkF（zk）‖≤2‖dk‖+αk‖d︿k‖+L1‖dk+αkd︿k‖2‖（JTkJk+λkI）－1JTk‖≤c3dist（xk，X*）

引理3 在假設1成立下，若xk∈N（x，05b），有[10]

（a）‖U1UΤ1Fk‖≤O（‖xk－x-k‖）；（b）‖U2UΤ2Fk‖≤O（‖xk－x-k‖2）；（c）‖U3UΤ3Fk‖≤O（‖xk－x-k‖2）

引理4 在假設1成立下，若xk，yk∈N（x，05b），有[10]

（a）‖U1UT1F（yk）‖≤O（‖xk－x-k‖2）；（b）‖U2UT2F（yk）‖≤O（‖xk－x-k‖3）；（c）‖U3UT3F（yk）‖≤O（‖xk－x-k‖3）

引理5 在假設1成立下，若xk，yk，zk∈N（x，05b），有[13]

（a）‖U1UΤ1F（zk）‖≤O（‖xk－x-k‖3）；（b）‖U2UΤ2F（zk）‖≤O（‖xk－x-k‖4）；（c）‖U3UΤ3F（zk）‖≤O（‖xk－x-k‖4）

定理1 當假設1成立，AMLM算法產生的xk四階收斂于式（1）的解。

證明：當δ∈［1， 2］，由引理1（a），引理5，式（11）及文獻[13]中式（29），（76），（79），得

‖d～k‖≤O（‖x-k－xk‖3）+O（‖x-k－xk‖5－δ）≤O（‖x-k－xk‖3）

‖F（zk）+βkJkd～k‖≤O（‖x-k－xk‖3+δ）+O（‖x-k－xk‖4）+O（‖x-k－xk‖4）≤O（‖x-k－xk‖4）

再結合式（4），式（5），式（6），式（7），及引理3，引理4及引理5，有

c‖x-k+1－xk+1‖≤O（‖xk－x-k‖4）（14）

此外，‖xk－x-k‖=dist（xk，X*）≤ ‖x-k+1－xk‖≤ ‖x-k+1－xk+1‖+‖sk‖，結合式（14），有‖x-k－xk‖≤2‖sk‖對充分大的k成立，再根據引理3，引理4及引理5，得‖dk+1+αk+1d︿k+1+βk+1d～k+1‖≤O（‖dk+αkd︿k+βkd～k‖4），即當δ∈［1， 2］，‖sk+1‖≤O（‖sk‖4），綜上，AMLM算法產生的xk四階收斂于式（1）的解。

3 數值實驗

實驗時，AMLM算法參數設置：p0=10－4，p1=025，p2=075，N0=5，m=10－8，μ0=10－4，η=10－4，δ=2，停機準則為‖JΤkFk‖≤10－5或迭代次數超過100（n+1），并與文獻[13]中提出的算法（記NMLM），文獻[9]提出的算法（記ALLM）以及文獻[17]提出的算法（記NALM）進行對比。NMLM算法參數：p0=10－4，p1=025，p2=075，m=10－8，μ0=10－4，δ=2；ALLM算法參數：p0=10－4，p1=025，p2=075，N0=5，m=10－8，μ0=10－4，θ=05，δ=2；NALM算法參數：q0=10－4，q1=025，q2=075，N0=5，m=10－8，μ0=10－4，α～=5，δ=1，ρ=05；具體結果見表1，其中x0取自文獻[19]，“n”表示函數的維數，“NF”表示函數的計算次數，“NJ”表示雅可比矩陣的計算次數，“NT=NF+NJ*n”表示總的計算次數，第四列表示初始點為－10x0，－x0，x0，10x0，100x0。實驗是在軟件平臺MatlabR2022b且配置為i7-7500U CPU，64位270GHz的個人計算機上執行。

測試問題根據文獻[20]中的方法修改而來，具體形式為

F︿（x）=F（x）－J（x*）A（ATA）－1AT（x－x*）

其中，F（x）是文獻[19]中給出的非奇異測試函數，x*是F（x）=0的解，A∈Rn×k（1≤k≤n）是列滿秩矩陣。顯然，F︿（x*）=0，J︿（x*）=J（x*）（I－A（ATA）－1AT）的秩為n－k，取A∈Rn×1，AT=（1，1，…，1），那么J︿（x*）的秩為n－1。

由表1可知，對于絕大部分測試問題的實驗結果，ALLM算法的函數的計算次數雖然較少，但是雅可比矩陣的計算次數和總的計算次數較多，而AMLM算法的雅可比矩陣的計算次數和總的計算次數相對較少，均少于NALM算法和ALLM算法；當算例2取初始點為x0、10x0、100x0，算例5取初始點為x0、100x0，算例4、算例6、算例7、算例9取初始點為10x0、100x0，算例3、算例8取初始點為100x0時，AMLM算法的函數的計算次數、雅可比矩陣計算次數、總的計算次數均小于NMLM算法。

4 結論

本文結合非單調技術及信賴域技術提出了一種求解非線性方程組的有效的LM算法（AMLM算法），在局部誤差界條件下，證明了AMLM算法具有局部快速收斂性。數值實驗結果表明，AMLM算法的結果具有一定的提升，并且穩定有效。然而在求解更大規模的問題時，如何提升收斂速度和節省更多的計算量在今后有待解決。

參考文獻

[1]LEONOV E A，POLBIN A V. Numerical search for a global solution in a two-mode economy model with an exhaustible resource of hydrocarbons[J]. Mathematical Models an Computer Simulations，2022，14（2）： 213-223.

[2]XU D，BAI Z Y，JIN X L，et al. A mean-variance portfolio optimization approach for high-renewable energy hub[J]. Applied Energy，2022，325：119888.

[3]MANZOOR Z，IQBAL M S，HUSSAIN S，et al. A study of propagation of the ultra-short femtosecond pulses in an optical fiber by using the extended generalized Riccati equation mapping method[J]. Optical and Quantum Electronics，2023，55（8）：717.

[4]YAMASHITA N，FUKUSHIMA M. On the rate of convergence of the Levenberg-Marquardt method[J]. Computing， 2001，15：239-249.

[5]FAN J Y. A modified Levenberg-Marquardt algorithm for singular system of nonlinear equations[J]. Journal of Computational Mathematics，2003，21（5）：625-636.

[6]FAN J Y，YUAN Y X. On the quadratic convergence of the Levenberg-Marquardt method without nonsingularity assumption[J]. Computing，2005，74：23-39.

[7]MA C F，JIANG L H. Some research on Levenberg-Marquardt method for the nonlinear equations[J]. Applied Mathematics and Computation，2007，184（2）：1032-1040.

[8]ZHENG L，CHEN L，MA Y F. A variant of the Levenberg-Marquardt method with adaptive parameters for systems of nonlinear equations[J]. AIMS Math，2022，7（1）：1241-1256.

[9]韓揚，芮紹平. 一種求解非線性方程組的修正Levenberg-Marquardt算法[J]. 青島大學學報（自然科學版），2023，36（1）：8-14.

[10] FAN J Y. The modified Levenberg-Marquardt method for nonlinear equations with cubic convergence[J]. Mathematics of Computation，2012，81（277）：447-466.

[11] FAN J Y. Accelerating the modified Levenberg-Marquardt method for nonlinear equations[J]. Mathematics of Computation，2014，83（287）：1173-1187.

[12] YANG X. A higher-order Levenberg-Marquardt method for nonlinear equations[J]. Applied Mathematics and Computation，2013，219（22）：10682-10694.

[13] CHEN L. A modified Levenberg-Marquardt method with line search for nonlinear equations[J]. Computational Optimization and Applications，2016，65（3）：753-779.

[14] AHOOKHOSH M，AMINI K. A nonmonotone trust region method with adaptive radius for unconstrained optimization problems[J]. Computers & Mathematics with Applications，2010，60（3）：411-422.

[15] AHOOKHOSH M，AMINI K. An efficient nonmonotone trust-region method for unconstrained optimization[J]. Numerical Algorithms，2012，59（4）：523-540.

[16] AMINI K，ROSTAMI F，CARISYI G. An efficient Levenberg-Marquardt method with a new LM parameter for systems of nonlinear equations[J]. Optimization，2018，67（5）：637-650.

[17] REZAEIPARSA Z，ASHRAFI A. A new adaptive Levenberg-Marquardt parameter with a nonmonotone and trust region strategies for the system of nonlinear equations[J/OL]. [2023-11-03]. Mathematical Sciences，2023. https：//webofscience.clarivate.cn/wos/alldb/full-record/WOS：000976510800001.

[18] STEWART G W，SUN J G. Matrix perturbation theory[M]. Boston：Academic Press，1990.

[19] MORE J J，GARBOW B S，HILLSTROM K E. Testing unconstrained optimization software[J]. ACM Transactions on Mathematical Software （TOMS），1981，7（1）：17-41.

[20] SCHNABEL R B，FRANK P D. Tensor methods for nonlinear equations[J]. SIAM Journal on Numerical Analysis，1984，21（5）：815-843.

An Effective Levenberg-Marquardt Algorithm for Solving Systems of Nonlinear Equations

HAN Yang，RUI Shao-ping

（School of Mathematical Sciences，Huaibei Normal University，Huaibei 235000，China）

Abstract：

A new effective Levenberg-Marquardt algorithm for solving systems of nonlinear equations，namely AMLM algorithm was presented by modifying LM parameters， combining nonmonotone technique and trust region technique. The local fast convergence of the AMLM algorithm was proved under the local error bound condition. Numerical results show that the AMLM algorithm is stable and effective.

Keywords：

LM parameter；nonmonotone technique；trust region technique；Levenberg-Marquardt algorithm