《非线性无约束优化问题的伪轨线追踪算法.doc》由会员分享，可在线阅读，更多相关《非线性无约束优化问题的伪轨线追踪算法.doc（11页珍藏版）》请在三一文库上搜索。
1、精品论文非线性无约束优化问题的伪轨线追踪算法刘舒天 北京邮电大学信息工程学院，北京(100876) E-mail: 摘要：无约束优化问题与常微分方程组的求解有着紧密的联系。本文分析了一种基于信赖域理论的伪轨线追踪算法，这类算法也可以被视为是具有自适应步长调节能力的线性化隐式 欧拉法。该算法的设计涵盖两个重要的数值分析领域，即最优化问题与常微分方程的数值解法。从最优化理论的角度来说，算法中使用了信赖域方法，同时，它受到了 Levenberg-Marquardt 参数设置的启发。这样可以使算法避免在每一步迭代中都必须满足信赖 域半径的约束。另一方面，从常微分方程的数值方法角度来说，算法构造的基础是
2、线性化隐 式欧拉法，这是算法加速的关键。本文中，我们分析了该算法的全局收敛性。进一步，我们 将算法应用于求解非线性无约束优化问题。数值算例显示，伪轨线追踪算法相较于传统的常 微分方程算法，在稳态点求解方面有着明显的优势。关键词：非约束优化；常微分方程组；信赖域；全局收敛性 中图分类号：01-8，0241.引 言本文的研究内容包含了两个重要的数值分析领域，即最优化问题与常微分方程的数值解 法。我们首先说明这二者之间的紧密联系。给定一个平滑函数 f ，一个有效的非约束优化算法旨在寻找一个局部最小点 x* ，使得对 于 x* 邻区内的任意点 x ，都满足 f ( x* ) f ( x) 。这等价于我
3、们考虑求解如下的非约束优化问题- 4 -min f ( x)xR n(1.1)针对式（1.1）所描述的非约束问题，我们有许多方法1, 2, 6, 9可以将其转化为如下的具有初始 值的常微分方程组x (t) = F ( x(t),t 0,x(0) = x0(1.2)其中，函数 F 连续可微，且 F ( x) = f ( x) 。下面给出的定理，指出了 x* 成为局部最小点的充分必要条件。定理的证明，可以在相关的书籍中4查找。定理 1.1 假设 x* 为局部最小点，那么其成立必要条件为f ( x* ) = 0 ，且 2 f ( x * ) 为半正 定的。进一步， x* 为局部最小点的充要条件，为
4、2 f ( x * ) 为正定的，且f ( x* ) = 0 。另一方面，根据高等数学中的“链式规则”，我们可以得到如下的表达式mdf ( x(t) dt = (f=1xi ) (ximit) = (fi =1x ) 2 = f ( x(t) 2 0(1.3)i显然，由式（1.3）我们可以发现，沿着该常微分方程组的任意解前进，随着时间t 的演进， 函数 f ( x(t ) 在二范数意义上总是不增的。进一步，若 f ( x(t) 0 ，则是严格递减的。所以，寻找具有初始值的常微分方程组（1.2）在t 足够时的解，可以被视为寻找非约束优化问题（1.1） 的局部最小值。本文从常微分方程的数值方法出发
5、，结合最优化理论中的方法，提出了一种新的算 法“伪轨线追踪算法”，用于求解非线性无约束优化问题。对于该算法的分析，我们可以 从两个不同的角度来思考。一方面，从最优化理论的角度来说，算法中使用了信赖域方法， 同时，它受到了 Levenberg-Marquardt 参数设置的启发。这样的处理手段，得以使算法避免 在每一步迭代中都必须满足信赖域半径的约束。另一方面，从常微分方程的数值方法角度来说，首先，算法构造的基础是隐式欧拉的迭代形式，进一步，我们又将其线性化。这是算法加速的关键。 接下来，本文将按照如下的顺序进行分析陈述。在下一节中，伪轨线追踪算法的具体实现方法以及算法分析将会被详细陈述。在第三
6、节中，我们将对算法的全局收敛性进行分析。 在第四节实际的数值算例中，我们将对比伪轨线追踪算法与传统的常微分方程算法的性能差 异。最后，是对全文进行的总结。2.伪轨线追踪算法2.1 伪轨线追踪算法的背景为了寻找函数 f 的局部最小值点，大多数的数值方法都是通过从一个初始的猜测值 x0 开始，产生一个搜索序列xk 。与之类似的，求解常微分方程组（1.2）的一步法同样产生一个 搜索序列xk 且 xk x(tk ) 。而时间点tk 则是动态的由时间步长所决定，即 tk = tk +1 tk 。一方面，对于经典的最速下降法，我们并不陌生。它具有如下的形式：xk +1 = xk k f ( xk )(2.
7、1)其中，k是一个常量，通常它可以通过线搜索得到。事实上，式（2.1）等价于将显式欧拉法应用于具有梯度形式的常微分方程组，且其时间步长 tk k。我们很清楚，当函数存在“陡峭山谷”，最速下降法的表现不尽如人意。这可以与显式欧拉法处理刚性问题的情形相类比。 另一方面，最优化理论中的牛顿法，是基于如下的局部二次模型：qk ( ) = f ( xk) + f ( xk)T + 1 T 2 f ( x )2k运用泰勒展式的相关知识，我们可以发现， qk ( ) 实际上是 f ( xk + ) 的二次逼近。所以，如果2 f ( x ) 是正定的，那么 q ( ) 具有唯一的最小值点k( 2 ()1( )
8、 。由此，我们可k k以得到牛顿法的表达式： k = f xkf xkkxk +1= xk ( 2 f ( x) 1 f ( x )(2.2)伪轨线追踪算法，正是基于上面这两个方面的研究。首先，我们将隐式欧拉法应用于常k微分方程组（1.2），且时间步长为tk，得到如下表达式xk +1 = xk t k f ( xk +1 )显然，这是一个非线性方程组。进一步，我们将牛顿法运用于式（2.3），得到2(2.3)整理得到xk +1 = xk t k f ( xk ) t k f ( xk )(xk +1 xk )(2.4)k其中 xk = xk +1 xk 。(t 2 f ( x) + I )xk=
9、 t k f ( xk )(2.5)通常，这一方法被称作线性化的隐式欧拉法3, 5。它具有非常好的性质。一方面，当很大时，式（2.5）变为tkkkxk +1 xk ( 2 f ( x) 1 f ( x )这是一个类似牛顿法（2.2）的表达式。另一方面，当xk +1 xk t k f ( xk )tk趋向于 0 时，式（2.5）蜕变为这恰好是具有小步长的最速下降法（2.1）。所以，（2.5）这一常微分方程的数值算法，在极限的情况下与经典最优化算法表现相当。然而，如果使用经典的局部误差原则去控制时间步长tk，方法（2.5）需要消耗大量不必要的计算。为了避免这种情况的发生，同时又精确地搜索到常微分方
10、程组（1.2）的稳态 解，我们使用信赖域的方法去控制步长，以期得到迫近平衡点时的快速收敛性。综上，将线性化的隐式欧拉法与信赖域的方法结合，我们提出如下的伪轨线追踪算法。在我们的算法中，使用 hk替代1 tk ，从而得到如下表达式k( 2 f ( x) + hk I )xk= f ( xk )(2.6)2.2 伪轨线追踪算法的实现方法算法 1 伪轨线追踪算法（PTM）步骤 0 参数初始化：选择初始猜测点，初始参数，以及常量 、 、 和 ，且满足设置 k = 00 1x0 1 1 和2 20 1h0 1 21 2 1 22步骤 1 定义逼近模型：计算 fk = f ( xk ) ，gk = f (
11、 xk ) ，Gk = f ( xk ) ，并定义逼近函数为k k kq (s) = f + sT g + 12ksT G s步骤 2 计算搜索步长：如果 Gk + hk I 正定，求解(Gk + hk I )sk = g k得步长 sk 与测验点 xtrial= xk+ sk，否则，设置 k= 1，去步骤 4。步骤 3 计算比率：设置 k =f ( xk ) f ( xtrial )f ( xk ) qk (sk )步骤 4 接受测验点：如果 k 0 ，那么设置 xk +1 = xk 。否则，设置 xk +1 = xtrial 。步骤 5 调整参数 hk：设置10hk , if 2 hk ,
12、 if k 00 k 1设置 k = k + 1 ，去步骤 1。hk = hk ,if 1 hk ,if1 0(3.2)由于 sk是局部子问题（2.8）的解，从而有qk ( sk s) qk (sk )。所以，qk = f ( xk ) qk (sk ) qk (0) qk ( sk s) 。进一步，Tqk sk sg k + o( sk ) = sk + o( sk ) 。(3.3)与此同时，根据泰勒展式2f ( xk ) = qk + o( sk )(3.4)根据式（3.1）、（3.3）、和（3.4），我们可以得到当 k , k S ， k = 1 + o(1) 。这与 k 1相矛盾。所以
13、， g = 0 。假设G = G( x ) 是非半正定的，那么存在一个归一化的方向 v ，使得v T G v = , 0(3.5) 取 k S 以及 = 1 ，使得v Tg k 0, k S , k k 。由于 sk 是局部子问题（2.8）的解，从而所以，212qk qk (0) qk ( sk v) sk2Tv Gk v(3.6)kq 1 sk2 k2 + o( s )(3.7)根据式（3.1）、（3.4）、和（3.7），我们可以得到当 k , k S ， k = 1 + o(1) 。这与 k 1相矛盾。所以G 是半正定的。情形（2）：根据式（2.7）中 hk的形式，可以推断必然存在一个无限
14、子序列，且它的索 引形成集合 S ，且 S S , k 2 , min (Gk + hk I ) 。假设 g = f ( x ) 0 ，那么，存在常量 gmin 0 和足够大的 k S ，使得 g k g min。所以g min g k Gk + hk I sk (Gmax + W ) sk其中，Gmax= supxB 2 f ( x) 。从而可以得到sk g min，足够大的 k S 。所以，Gmax + WGinf l k lmin =+ Wg min(3.8)k S maxk令 f = f ( x ) ， f= f ( xk) f ( x1k +1) 0 。由 f f fkS，可以得到当
15、kk , k S 时，f k 0 。由于 k 2 ，所以 qk 0 。令 q(s) = f + s T g + 1 2 s TG s ，取 l (0, l min ) ，令 s 使得 q (s) 在 s h 上最小化，同时设 x = x + s 。那么，对于足够大的 k S ，有x xk s +x x= s + o(1) l + o(1) lk所以，qk ( x xk ) qk (sk ) = f k qk从而， g 0 是矛盾的。进一步，由于 g = 0 ，我们得到(3.9)ksk (Gk+ h I ) 1g 1 gk k 0 ，足够大的 k Sk假设G 不是半正定的。那么，式（3.5）与（
16、3.7）成立，且 = 1 + o(1) 。从而，由式（2.7）得 hk 0 。又因为 min (Gk + hk I ) ，必然有G 是半正定，矛盾产生。证毕4.数值实验我们通过两个特征值相关的数值实验，用伪轨线追踪算法（PTM）与 Matlab 常微分方 程组算法库中方法作对比。通过 Rayleigh 商的定义，我们可以将特征值问题转化为最优化 问题。如下：minxR nx T Ax x T x通过此最优化问题，我们可以寻找到矩阵 A 的最大与最小特征值。实验环境为一台装有 Matlab 的双精度浮点运算的 PC，中央处理器为 Intel(R) Core(TM)2CPU 4400 2.00GH
17、z。所有的 ODE 算法均在相同条件下终止，即F ( xk ) 其中， 是一个预设值。在我们的实验中，我们设为 1e-6。同时，我们为伪轨线追踪算法（PTM）选取1 = 0.25 ， 2 = 0.8 ， 1 = 0.5 ， 2 = 2 ，以及 h0 = min f ( x0 ) ,100 ，为 Matlab 常微分方程组算法设置 RelTol=1e-6 和 AbsTol=1e-9。有关常微分方程组的解法以及 Matlab 常微分方程组算法库的信息，可以查阅参考文献7, 8。4.1 实验一实验一是一个稠密矩阵，它包含有三个相隔的特征值簇。我们可以按照如下的方式构建：1) 选择 = diag(1,
18、1,0,0,1,.,1) R nn2) 令C = rand (n, n) ，同时Q, R = qr(C)3) 定义 A = Q T Q实验中，我们选取 n = 1000 。Extreme eigenvalue problem0-1PTM ODE45-2log|Ax - x|10-3-4-5-60 1020 304050 60708090number k of iterations图 1 具有分离的特征值簇的稠密矩阵由图一，我们可以发现，伪轨线追踪算法（PTM）比 ODE45 的收敛速度快几乎两倍。与此同时，PTM 具备更高的精度。由表一的数据里，同样可以得到这样的结论。表 1 实验一计算量算法
19、PTMODE45CPU 时间（s）1.717.4迭代次数41874.2 实验二实验二与实验一类似，不同之处在于两个特征值簇紧挨在一起。所以，需要算法有很高 的灵敏度。我们用实验一的方式构建，并将第一步替换为1) 选择 = diag(1e 4,1e 4,0,0,1,.,1) R nn实验中，我们同样选取 n=1000。从图二中，我们发现了与实验一中类似的情形。然而，在此实验中，ODE45 遇到了一 个差错平底。PTM 则克服了这一缺陷，这是由于 PTM 使用了信赖域的方法来控制步长，而 ODE45 采用的传统的差错控制理论。表二用数据的形式反映了图二的现象。Extreme eigenvalue
20、problem0-1-2-3log|Ax - x|10-4-5PTM ODE45-70 50 100 150200250300number k of iterations图 2 具有紧靠着的特征值簇的稠密矩阵表 2 实验二计算量算法PTMODE45CPU 时间（s）2.2停滞迭代次数51停滞5.结论本文给出了一种新的用于求解非线性无约束优化问题的算法 伪轨线追踪算法（PTM）。它是线性化的隐式欧拉法与信赖域方法的有机结合。一方面，它避免了去复杂计 算一个非线性方程组，另一方面，通过信赖域方法，算法本身得以加速。数值算例显示，该 算法在处理稳态解问题时，优于传统的常微分方程数值解法。原因在于，该
21、算法不考虑中间 态，而是通过伪轨线加速追踪最终的稳态，以获得速度的提升，又不失精确度。参考文献1 C. A. Botsaris. A class of methods for unconstrained minimization based on stable numerical integration techniques J. Journal of Mathematical Analysis and Applications, Vol. 63, pp. 729-749, 1978.2 A. A. Brown and M. C. Bartholomew-Biggs. Some effecti
22、ve methods for unconstrained optimization based on the solution of systems of ordinary differential equations J. Journal of Optimization and Theory applications. Vol.62, pp. 211-224, 1989.3 P. Deuflhard. Newton Methods for Nonlinear Problems: Affine Invariance and Adaptive Algorithms M. Spring-Verla
23、g, 2004.4 R. Fletcher. Practical Methods of Optimization M. 2nd ed., John Wiley,Chichester, 1987.5 E. Hairer and G. Wanner. Solving Ordinary Differential Equations II, Stiff and Differential-Algebraic ProblemsM. 2nd ed., Springer-Verlag, Berlin, 1996.6 J. Schropp. A note on minimization problems and
24、 multistep methods J. Numer. Math., Vol. 78, pp. 87-101,1997.7 L. F. Shampine and M. W. Reichelt. The Matlab ODE Suite M. SIAM J. Sci. Comput., Vol. 18, pp. 1-22,1997.8 刘德贵,费景高,韩天敏.刚性大系统数字仿真方法M. 河南科学技术出版社,1996.1.9 袁亚湘,孙文瑜.最优化理论与方法M.科学出版社,2003.4.Pseudo- trajectory method for nonlinear unconstrained o
25、ptimization problemLiu ShutianSchool of Information Engineering, Beijing University of Posts and Telecommunications, Beijing(100876)AbstractUnconstrained optimization problems are closed related to ordinary differential equations with gradient structure. In this paper, one trust region type method c
26、alled pseudo-trajectory algorithm is developed.This kind of algorithm can be regarded as linearized implicit Euler method mixed with adaptive trust region time step control. From the view of optimization theory, the algorithm utilizes trust-regionmethod. Furthermore, it is inspired by the configurat
27、ion of Levenberg-Marquardt parameter, which avoids satisfying the constraint of trust-region radius at each step. From the view of ordinarydifferential equation, the algorithm is based on the linearized implicit Euler method which leads to the acceleration of our method. In the work, we also provide
28、 the global convergence analysis. Furthermore,we apply this algorithm to solve nonlinear unconstrained optimization problem. Numerical experiments indicate that the new method is applicable in the problem of finding equilibrium steady state.Keywords: unconstrained optimization; ordinary differential equation; trust region; global convergence