一种求解非线性方程组的有效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