一类修正的阻尼牛顿法

一类修正的阻尼牛顿法

庞军彦, 李秦

（兰州交通大学 数理与软件工程学院，兰州 730070）

摘要：本文在Marquardt –Levenber 方法和Goldstein –Price 方法的基础上对阻尼牛顿法(1)()2()1()[()]()k k k k k x x f x f x λ+-=-??作了适当改进，得出了一种新的算法，与原来算法相比较，新算法避免了二阶导数矩阵的奇异性和非正定性，从而使迭代在二阶导数矩阵奇异和非正定的条件下也能进行．文章还给出了新算法收敛性分析和算法步骤，最后给出了数值试验．

关键字：阻尼牛顿法；迭代格式；收敛性；数值试验

A modified damped Newton method

PANG Jun-yan ，LI Qin

(School of Mathematics, Physics and Software Engineering, Lanzhou Jiaotong University,

Lanzhou 1730070, China)

Abstract: In this paper,based on the two methods of Marquardt –Levenber and Goldstein –Price.The damped Newton method (1)()2()1()[()]()k k k k k x x f x f x λ+-=-?? is improved properly.We proposed a new algorithm.In comparison with the original algorithm.The new algorithm avoids the singular and non positive definiteness of two order derivative matrix.The iteration can continue in the singular and non positive definite conditions of two order derivative matrix. The convergence of the new algorithm is analyzed and the algorithm steps is given.Finally,the numerical tests are presented.

Key words: damped Newton method; iterative format; convergence; numerical tests

1、 引言

在科学和工程计算中，常常要用到数值计算方面的知识来求最优化极小值

基金项目：国家自然科学基金（A020205)

作者简介：庞军彦（1983-），男，甘肃省陇南市礼县人，兰州交通大学硕士研究生，主要从事数值计算方

面的研究

李秦（1965-），女，陕西西安人，教授，主要研究方向为算法设计与分析

问题min ()f x ，而在一系列数值计算方法中，牛顿迭代法无疑是求此类问题的一种非常重要的方法，其迭代格式为(1)()2()1()[()]()k k k k x x f x f x +-=-?? （k=1,2，…） 其中()()k f x ?是函数()f x 在()k x 的一阶导数矩阵，2()()k f x ?是函数()f x 的二阶导数矩阵（Hessian 矩阵），众所周知，牛顿法有一定的缺点，即对初值要求很苛刻，要求初值0x 必须在*x 附近取值，否则有可能迭代不收敛或收敛很慢 ．阻尼牛顿

法是对牛顿法的修正得到的，其迭代格式为(1)()2()1()[()]()k k k k k x x f x f x λ+-=-??，其中k λ满足()()()()min k k k k k f x p f x p λ

λλ+=+． 阻尼牛顿法保持了牛顿法的二阶收敛性，即满足(1)*()*k k x x M x x +-≤-，而且对初值没有过于苛刻的要求。本文就是通过对阻尼牛顿法的修正，得出了一种新算法，并进行了收敛性分析．

2、新算法的得出

阻尼牛顿法的迭代格式为(1)()2()1()[()]()k k k k k x x f x f x λ+-=-?? （1） 其中k λ满足 ()()()()min k k k k k f x p f x p λ

λλ+=+．

在（1）中要求算 2()1[()]k f x -?，因而要保证二阶导数矩阵2()()k f x ? 必须是非奇异的，否则 2()1[()]k f x -?不存在，迭代就进行不下去，甚至有时候算2()1[()]k f x -?也很困难，在（1）中还要要求2()()k f x ?是正定的，否则迭代也不能够产生新点．鉴于以上缺点，本文在阻尼牛顿法[1],[2],[3]，Marquardt —Levenberg 法[1]和Goldstein —Price 法[1]三种方法的基础上对阻尼牛顿法进行修正．

在（1）中我们用()()k M x 近似代替二阶导数矩阵2()()k f x ?，则迭代格式变为 (1)()()1()[()]()k k k k k x x M x f x λ+-=-? （2） 其中k λ满足 ()()()()min k k k k k f x p f x p λ

λλ+=+，()1()[()]()k k k k p M x f x λ-=-?．

（2）就是本文新算法的迭代格式．新算法的关键是()()k M x 的得出， ()()()()k k M x Q x I =+?，其中I 为单位矩阵，?为一正常数，()()k Q x 为一矩阵， ()()()()12()((),(),......,())k k k k n Q x Q x Q x Q x =，()()k i Q x 求法如下：

()()()()()

()k k k k i i k f x e f x Q x ββ?+-?=，其中i e 是第i 个分量为1的单位向量，即

1,……i ……n

i e =(0,…..1……0), 而k β为：0βγ=，1,1,2,...k k p k βγ-==，γ为一给定常数．

在（2）中，由于()()()()k k M x Q x I α=+，故

① 当()()k Q x 正定时，取α=0 则()()()()k k M x Q x =，迭代格式就变为

(1)()()1()[()]()k k k k k x x Q x f x λ+-=-? （3） ② 当()()k Q x 不正定时，将?的值取得足够大，可认为是第二项I ?在起主导作用，从而能保证()()k M x 的正定性，这样既避免了原阻尼牛顿法中2()()k f x ?的非正定和奇异，而且还可避免二阶矩阵的计算．

3、 收敛性分析及算法步骤

3.1收敛性分析

定义1：设迭代过程(1)()()k k x x ?+=收敛于方程()x x ?=的根*x ，如果迭代误差*()k k e x x =-，当k →∞时成立1k p

k e c e +→(c 为常数且c ≠0)，则称迭代过程是p 阶收敛的。特别地，当p=1时称为线性收敛的，当1＜p ＜2时称为超线性收敛的，当p=2时称之为平方收敛的[2]．

定义2：设F 为从集合n D R ?到m R 单值映射，设x D ∈，若存在[0,)K R +∈=∞，满足''()(),p F x F y K x y -≤- ,x y X ?∈ ，(0,1]p ∈ ，K 为Lipschitz 常数，则称'F 是以K 为常数的Lipschitz 连续算子[6]．

引理1：设()f x 在点*x 可微，如果存在方向d ，使*()0T f x d ?<，则存在数0δ>，使得对每个(0,)λδ∈，有**()()f x d f x λ+<[3]．

引理2：设f 是n R 上的二阶连续可微函数，2()f x ?在点*x 的一领域内满足Lipschitz 条件（其常数是L ），且2()f x ?正定，其中*x 是满足二阶充分条件的局部极小值点，设{}k x 由经典牛顿法生成，则{}k x 二次收敛到*x [4]．

定理1：设()M x ：n n D R R ?→是Q —可导算子，且()M x 是2()()k f x ?的近似矩阵，()M x 满足p —Lipschitz 连续()()p

M x M y K x y -≤-，(0,1]p ∈，,n x y R ?∈，K 为正常数，且*1*()M x B -=,在阻尼牛顿算法中的阻尼常数

0k λε≥>,*x 为()F x ?的单零点，*(,)O x r D

?表示开球域，其中111

11().(1)2(*)p P

P p r p B K δ+=-+,(01)δ<<,则对任意*0(,)x O x r ∈，阻尼牛顿法(1)()()1()[()]()k k k k k x x M x f x λ+-=-?所产生的迭代序列{n x }有意义，且*n x x →． 证明：对定理条件知，对任意的,x y D ?∈有()()p

M x M y K x y -≤-，故更有*()()p M x M x K x y -≤-，取任一*(,)x O x r ∈，因而*x x r -< 故*1*12p

p B K x x p +-<<+且**1*0()1*p B M x B K x x -<--，则1x 有意义，且 **110000[()]()x x x x M x F x λ--=--?

=*1**0000000(1)()[()][()()()()]x x M x F x F x M x x x λλ---+?-?-- 从而**1**10000000(1)().()()()()x x x x M x F x F x M x x x λλ--≤--+?-?--

*1**0000*0(1).11*p p B K x x x x p B K x x λλ+≤--+-+-- =*0*000*0*1.(1*)(1)p p B K x x x x B K x x p λλ??- -+?- ?--+?

? 由于*0x x r -<，从而*01*(1)2p p B K x x p

λ+-<-+，所以 001(1)211(1(1))(1)2p p p p p

λλλλ+-+-++--++00111(1)p δλλδ-=-+++011[1]1(1)p δλδ-=+-++ 0(2)11(1)p p δλδ-+=+++α<<(2)11(1)p p δεδ

+-++<1 **10x x x x r αα-<-<

所以*1(,)x O x r ∈，由此类推*(,)n x O x r ∈，1n x +有意义且

***10...0n n n n x x x x x x r ααα+-<-<<-<→，从而*n x x →．

定理2：设F 是n R 上的二阶连续可微函数，()M x 是2()()k f x ?的近似矩阵，

若()M x 在点*x 的附近满足()()k M x M y L x y λ-≤-，,n x y R ∈，

1*1()()k M x M x β--≤（L 、β是常数）

，其中*x 是满足二阶条件的极小值点，{}k x 由（2）得到，则{}k x 二阶收敛到*x ．

证明：由于()M x 在点*x 附近满足()()k M x M y L x y λ-≤-，()M x 是2()()k f x ?的近似矩阵，所以**()()()()k k k k M x x x f x f x λ--?-?

2**()()()()k k k k f x x x f x f x λ≈?--?-?

122***0[()(())]()k k k k f x f x t x x x x dt λ=

?-?+--? 1

***0[()(())]()k k k k M x M x t x x x x dt λ≈-+--?

12*0(1)k L x x

t dt ≤--? 2*2

k L x x =- 所以得到1**()()k k k k k x x x M x f x x λ+-=-?-

1**()(()()()())

k k k k k M x M x x x f x f x λ-=--?+? 2*1*()2k L M x x x β-≤-

由定义1知{}k x 二阶收敛到*x ．

定理3：设函数()f x 在点*x 二阶连续可微，若*x 是()f x 的局部极小点，则*()0f x ?=．

证明：用反证法，设*()0f x ?≠，令方向1[()]()k d M x f x λ-=-?，其中()M x 是2()f x ?的近似矩阵，*()T f x d ?＝**1*()[()]()T k f x M x f x λ--??，由于()M x 正定且*

()0f x ?≠，

故**1*()[()]()0T f x M x f x -??>，所以*()0T f x d ?<，根据引理1，必存在0δ>，使得

当(0,)λδ∈时成立**()()f x d f x λ+<，这与*

x 是局部极小值点矛盾，所以*()0f x ?=．

3.2 算法步骤

给定(0)x ，正数r ，α，0ε>．

1、0k ?，r β?．

2、求()()k f x ?，()()k i f x e β?+，k r p β?，

()()()()()

()k k k k i i k f x e f x Q x ββ?+-?=，（i =1,2，…n ）．

3、记矩阵()()()()12()((),(),......,())k k k k n Q x Q x Q x Q x =．

4、若非奇异、正定，则做5，否则做6．

5、()1()[()]()k k k p Q x f x -=-?，转向7．

6、()1()[()]()k k k p M x f x -=-?,其中()()()()k k M x Q x I α=+（α为满足()()k M x 正定

的足够大的正常数，I 为单位矩阵）．

7、求(1)()k k k k x x p λ+=+，其中k λ满足()()()()min k k k k k f x p f x p λ

λλ+=+．

8、若(1)()k k x x ε+-≤，则求出了最优解*(1)k x x +=；若(1)()k k x x ε+->，则令

1k k +?转到2．

4、数值试验

例22121122(,)(1)f x x x x x x =+++

取()(0)0,0T x =，显然(0)0()2f x ???= ???，2(0)01()12f x ???= ???

2(0)1(0)

0[()]()p f x f x -=-??=20-?? ???，而(0)0()(2,0)f x p f λλ+=-216λ=的极小值点在0λ=达到，因此用阻尼牛顿法不能产生新点，迭代不能再进行下去，很明显，()(0)0,0T

x =不是12(,)f x x 的极小点，这是因为2()f x ?不是正定矩阵而引起的． 如果用本文得到的新算法进行计算：

(0)(0)(0)

0110()()()f x e f x Q x ββ?+-?==41?? ???（取01β=）

(0)(0)(0)

0220()()()f x e f x Q x ββ?+-?==12?? ???（取01β=） (0)(0)(0)12()((),())Q x Q x Q x =4112??= ???

正定，所以 (0)1(0)0[()]()p Q x f x -=-?48,77T

??=-- ???，而 (0)048()(,)77f x p f λλλ+=--4224328()(1)7497

λλλ=-++- 令(0)0()()f x p f λλ+=，从而求min ()f λ 令5

'

342419216()777f λλλ=+-=0，解之得0.56λ=，从而可以产生新点，迭代能继续进行．

参考文献

[1] 席少霖,赵凤治,最优化计算方法[M].社会科学技术出版社,1983.

[2] 孙凤芝.数值计算方法与实验[M].黑龙江大学出版社,2013.

[3] 陈开周.最优化计算方法[M].西北电讯工程学院出版社,1985.

[4] 吴佳,张艺.变分分析与优化[M].科学出版社,2013.

[5] 曹建胜,周武.牛顿法及带阻尼牛顿法的收敛域定理[J].南京师范大学学报(自然科学版),1989,12(2):24-27.

[6] 谢莹,张奠成.关于阻尼牛顿法收敛域的一个定理[J].高等学校计算数学学报,1986,1:85-87.

[7] 张建军,李春泉,张烈辉,影响阻尼牛顿法收敛性的两个重要参数[J].纯粹数学与应用数学,2012,28(4):433-439.

[8] 吕勇,刘兴国.牛顿迭代法加速收敛的一种迭代格式[J].武汉科技学院学报,2006,19(2):68-70.

[9] 黄争鸣.一种改进的牛顿法[J].1986,14(6):865-870.