非线性方程的不动点迭代方法研究

重庆文理学院

2011-2012下学期《数值方法》课程论文

题目：非线性方程的不动点迭代方法研究

中国﹒重庆

2012 年 06 月

通过从现实的一个球体的实际问题引出对非线性方程的不动点迭代研究，在理解迭代规则1()k k p g p +=的基础上通过对迭代法和不动点迭代法的基本思想即找()0f x =的同解变形()x g x =，然后运用初值0x 迭代，求出误差范围内的近似解1k x p +≈。运用函数连续性证明不动点的存在性和运用中值定理和均值定理证明不动点唯一性，进而推导出不动点迭代法的推导步骤。然后又用均值定理和数学归纳法证明出收敛性，并在此基础上引出误差边界。再通过对开题提出的球体问题案例的求解，进一步来加深非线性方程对不动点迭代法实证说明，由此联系到不动点迭代法在其他一些领域如物理和工程等的运用。

不动点迭代式（2.2）通常只有线性收敛，有时甚至不收敛，进而在原有的基础上拓展到加速迭代法的收敛性的讨论通常，从而对Steffensen 加速迭代和Aitken （埃特金）加速迭代的讨论。

关键字 不动点迭代； 收敛性； Steffensen 加速迭代； Aitken 加速迭代

1 问题的提出--------------------------------------------------------1

2 算法的思想--------------------------------------------------------1 2.1迭代法的基本思想----------------------------------------------1

2.2 不动点迭代法的基本思想----------------------------------------2

3 算法的推导及步骤--------------------------------------------------2 3.1 算法的推导----------------------------------------------------2

3.2 算法的步骤----------------------------------------------------3

4 算法的分析--------------------------------------------------------4 4.1 收敛性分析----------------------------------------------------4 4.2 误差性分析----------------------------------------------------6

4.3 稳定性分析----------------------------------------------------7

5 算法的实现--------------------------------------------------------7 5.1 案例----------------------------------------------------------7 5.2 求解过程------------------------------------------------------7

5.3 不动点迭代法代码及输出结果------------------------------------8

6 运用举例---------------------------------------------------------10

7 知识拓展---------------------------------------------------------10 7.1 Steffensen加速迭代-------------------------------------------10 7.2 Aitken（埃特金）加速迭代法-------------------------------------11

1 问题的提出

在现实生活当中我们会遇到很多关于诸如像球体的物理和工程问题，例如：

球体的半径为r ，并浸入水中，深度为d ，假设这个球由由一种密度为ρ=0.638的长叶松构成，且它的半径r =10cm 。当球浸入水中时，它的进水的质量为多少？而这些现实中的问题的解决都要涉及到求解方程的根、线性和非线性方程组以及微分方程的数值解。 如上例问题的解决：

当一个球以深度d 浸入水中时，所排开水的质量w M 为：

2

220(3)(())3

d x w d r d r x r d M ππ---==? （1.1）

而球体的质量为：

34/3b r M πρ= （1.2）

根据阿基米德定律有：

b w M M = （1.3） 由方程（1.1）、（1.2）、（1.3）联立得： 323(34)03

d d r r πρ-+=

而当r=10, ρ=0.638时，方程为：

23(255230

)

03

d d π-+=

此时这个物理学中的实际问题就变为了一个3次非线性方程23255230y d d =-+的根，在根据具体现实情况舍去不合理的根。

要解决如上的一些实际问题，就需要求解方程的根、线性和非线性方程组的解以及微分方程的解。计算机科学中的一个基本要素就是迭代，迭代技术用来求解方程的根、线性和非线性方程组的解以及微分方程的解，而不动点迭代法就是迭代法的一种基本方法。在此我们仅对非线性方程的不动点迭代进行研究。

2 算法的基本思想

2.1 迭代法的基本思想

迭代算法是数值计算方法中一种逐次逼近的方法，首先将方程()0f x =改写成某种等价的形式，由等价形式构造相应的迭代公式然后选取方程的某个初值近似

根o x 代入公式反复校正根的近似值，直到满足精度为止。

2.2 不动点迭法代的基本思想

将非线性方程()0f x =化为等价形式：

()x x ?= (2.1)

***()()0f x x x ?=?=

称*x 为函数()x ?的一个不动点，给定初始近似值o x ，可以得到0()x η?=，如此反复，构造迭代公式：

(1)(),0,1,2,3......k k x x k ?+== （2.2） 此处()x ?为迭代函数。如果对任何初值[]0,x a b ∈，由（2.2）得到的序列{}k x 有极限

*lim k x x x →∞

=

则称迭代公式（2.2）收敛，且**()x x ?=是()x ?的不不动点，故（2.2）就是不动点迭代法

3 算法的推导及步骤

3.1 算法的推导

3.1.1 不动点的存在性推导及证明

设g 是一个连续函数，且{}0n n p ∞

=是由不动点迭代生成的序列。如果

l i m n x p p →∞

=，则p 是()g x 的不动点。

证明：

如果lim n x p p →∞

=，则1lim n x p p +→∞

=。根据这个结论，g 的连续性和1()

n n p g p +=存在如下关系：

1()(lim )lim ()lim n n n x x x g p g p g p p p +→∞

→∞

→∞

==== （3.1）

因此，p 是()g x 的不动点。

3.1.2 不动点的唯一性推导及证明

设函数[],g C a b ∈。如果对于所有[],x a b ∈，映射()y g x =的范围满足

[],y a b ∈，则函数g 在[],a b 内有一个不动点。 （1）

此外，设'()g x 定义在(,)a b 内，且对于所有(,)x a b ∈，存在正常数1K <，使得，'1()g K x ≤<，则函数g 在[],a b 内有唯一的不动点p 。 （2） 证明：

对于命题（1）

如果()g a a =或()g b b =，则断言为真；否则()g a 必须满足()(,]g a a b ∈，()g b 的值必须满足()[,)g b a b ∈。表达式()()f x x g x ≡-有如下特性： ()()0f a a g a ≡-<且()()0f b b g b ≡->

对()f x 应用中值定理，而且由于常量0L =，可推断出存在数P ,且(,)P a b ∈，满足()0f P =。因此，()P g P =，且P 是'()g x 的不动点。

对于命题（2）

必须证明结果是唯一的。采用反证法，设存在两个不动点1P 和2P 。根据均值定理，可推断出存在数(,)d a b ∈，满足： '2121

()()()g P g P g d P P -=

- （3.2） 根据假设，有11()g P P ==且22()g P P ==，对3.2式的右边进行化简可得 '21

21

()1P P g d P P -=

=- 但是这与命题（2）中假设在(,)a b 内有'1()g x <矛盾，因此不可能存在两个不动点。所以，在命题（2）的假设条件下，()g x 在[],a b 内有一个唯一的不动点P 。

3.2 算法的步骤

（1）确定方程()0f x =的等价形式()x g x =，为确保迭代工程的收敛，要求

()g x 满足李普希茨条件（或'()1g x q ≤<）； （2）选取初始值0x ，按公式

(1)(),0,1,2,3......k k x g x k +==进行迭代； （3）若1k k x x ε+-<，则停止计算1k x x +≈。

4 算法分析

4.1收敛性分析

设有：

①g ，'[,]g C a b ∈； ②K 是一个正常数； ③0(,)p a b ∈；

④对于所有[,]x a b ∈，有()[,]g x a b ∈。

如果对于所有[,]x a b ∈，有'1()g K x ≤<，则迭代1()n n p g p -=将收敛到唯一的不动点[,]P a b ∈。在这种情况下，P 称为吸引（attractive ）不动点。 （3）

如果对于所有[,]x a b ∈，有'1()g x >，则迭代1()n n p g p -=将不会收敛到P 这种情况下，P 称为排斥（repelling ）不动点，而且迭代显示出局部发散性。（4）

批注：

①在命题（4）中假设0p P ≠。

②因为函数g 在包含P 的一段间隔中是连续的，可在命题（3）和命题（4）中分别利用更简单的判别条件'1()g K x ≤<和'1()g x >。 证明：

对于命题（3）

首先要证明点{}0n P ∞

=都位于(,)a b 内。从0p 开始，根据均值定理，可推导出存在一个值0(,)c a b ∈满足

1000'()()()()P p g P g p g c P p -=-=-

0000'()g P p K P p P p c =-≤-<- （3.3） 因此，1p 比0p 更接近P ，且1(,)p a b ∈如图1.一般情况下设1(,)n p a b -∈，则：

1011'()()()()n n P p g P g p g c P p ---=-=-

1111'()n n n n g P p K P p P p c ----=-≤-<- （3.4） 因此，(,)n p a b ∈，而且可归纳出所有的点{}0n n p ∞

=位于(,)a b 内。

图1、P ，0p ，1p ，0P p -和1P p -之间的关系

为完成命题（3）的证明，需要证明如下表达式成立：

lim 0n x P p →∞

-= （3.5）

首先，用归纳法的证明可建立如下不等式：

0n n P p K P p -≤- （3.6） 当1n =时满足关系式3.3。利用归纳假设10n n P p K P p --≤-和关系式3.4的思路，可得到：

110n n n n n P p K P p KK P p K P p ---≤-≤-=- 这样，通过归纳法可以得出，所有的n 满足不等式3.6.由于01K <<，所以当n 趋于无穷大时，项n K 趋近于0.因此

00lim lim 0n n n n P p K P p →∞

→∞

≤-≤-= （3.7）

0P p -的极限压缩在左右两边的0之间，所以可得出lim n n P P →∞

=，且根据前面定

理3.1，迭代1()n n p g p -=收敛到不动点P 。因此命题（3）得证。

对于命题（4）

首先要证明点{}0n P ∞

=都位于(,)a b 内。从0p 开始，根据均值定理，可推导出存在一个值0(,)c a b ∈满足：

1000'()()()()P p g P g p g c P p -=-=-

0000'()g P p H P p P p c =->-=-（3.8） 因此，0p 比1p 更接近P ，且1(,)p a b ∈，一般情况下设1(,)n p a b -∈，则：

1011'()()()()n n P p g P g p g c P p ---=-=-

1111'()n n n n g P p H P p P p c ----=->-=-（3.9） 因此，(,)n p a b ∈，而且可归纳出所有的点{}0n n p ∞

=位于(,)a b 内。为完成命题（3）的证明，需要证明如下表达式成立：

lim n x P p →∞

-=∞ （3.10）

首先，用归纳法的证明可建立如下不等式：

0n n P p H P p -≥- （3.11）

当1n =时满足关系式3.8。利用归纳假设10n n P p H P p --≥-和关系式3.9的思路，可得到：

110n n n n n P p H P p HH P p H P p ---≥-≥-=-

这样，通过归纳法可以得出，所有的n 满足不等式3.11.由于01K <<，所以当n 趋于无穷大时，项n K 趋近于∞。因此

0lim lim n n n n P p H P p →∞

→∞

-≥-=∞ （3.12）

0P p -的极限延伸到右边的∞之外，

所以可得出lim n n P →∞

=∞。因此命题（4）得证。

4.2 误差分析 设函数g 满足：

①g ，'[,]g C a b ∈； ②K 是一个正常数； ③0(,)p a b ∈；

④对于所有[,]x a b ∈，有()[,]g x a b ∈。

当用n p 去近似P 时，引入误差边界如下：

0n n P p K P p -≤- 对于所有1n ≥ 0

1n n K P p P p K

--≤- 对于所有1n ≥

4.3稳定性分析

在原有不动点迭代法的基础上将原有的初始值0x 修改为0x x +?然后按照迭代规则：

(1)(),0,1,2,3......k k x g x k +==

求出迭代通项最后比较其最后的所得值与真实值满足误差范围内内则其算法是稳定的即：

1k x x ε+-<

5 算法实现

5.1案例

在此案例就已论文开头提出的关于物理球体问题为例：球体的半径为r ，并浸入水中，深度为d ，假设这个球由由一种密度为ρ=0.638的长叶松构成，且它的半径r =10cm 。当球浸入水中时，它的进水的质量为多少？

5.2求解过程

当一个球以深度d 浸入水中时，所排开水的质量w M 为：

2

220(3)(())3

d x w d r d r x r d M ππ---==? （5.1）

而球体的质量为：

34/3b r M πρ= （5.2）

根据阿基米德定律有：

b w M M = （5.3） 由方程（5.1）、（5.2）、（5.3）联立得： 323(34)03

d d r r πρ-+=

而当r =10, ρ=0.638时，方程为：

23(255230

)

03

d d π-+=

此时这个物理学中的实际问题就变为了一个3次非线性方程23255230y d d =-+的根。对于此方程的求解可以用到不动点迭代法求解。解方程步骤如下： 首先有题意可知起[0,20]d ∈

①方程230255230d d =-+的等价形式为：

()d g d ==

当[0,20]d ∈时：

'

()1g d =

<

②选取初始值0x =12然后按公式(1)(),0,1,2,3......k k x g x k +==进行迭代。 ③若121(0.1)k k d d ε+-<=，则停止计算1k d d +≈。

5.3 不动点迭代法代码及输出结果 代码：

function [k,p,err,P]=fixpt(g,p0,tol,max1)

P(1)=p0;

for k=2:max1

P(k)=feval(g,P(k-1)); err=abs(P(k)-P(k-1));

relerr=err/(abs(P(k))+eps); p=P(k);

if (err

if k == max1

disp('maximun nuumber of iterations exceeded') end P=P';

[k,p,err,P]=fixpt(gx, 12, 0.000000000001, 200)

输出结果： k =

45

p =

11.86150150813511

err =

1.007904870675702e-011 P =

12.00000000000000

11.94431524477928

11.91085729014850

11.89085883104665

11.87894290969856

11.87185626654913

11.86764642816162

11.86514723509375

11.86366416540031

11.86278429075393

11.86226235240524

11.86195276626302

11.86176914523977

11.86166023954459

11.86159564865799

11.86155734083855

11.86153462122358

11.86152114671514

11.86151315529503

11.86150841577415

11.86150560487935

11.86150393780662

11.86150294910695

11.86150236273369

11.86150201497027

11.86150180872044

11.86150168639880

11.86150161385288

11.86150157082771

11.86150154531055

11.86150153017696 11.86150152120160 11.86150151587854 11.86150151272156 11.86150151084924 11.86150150973881 11.86150150908025 11.86150150868967 11.86150150845802 11.86150150832064 11.86150150823916 11.86150150819084 11.86150150816218 11.86150150814519 11.86150150813511

6 应用举例

不动点迭代法的应用极其广泛，在物理和工程学中都有相应的运用，如：不动点迭代法在齿形系数和应力修正系数计算中的应用

利用不动点迭代法求解危险截面齿厚和弯曲力臂为进一步用有限元法求解齿轮的弹性变形带来便利。

不动点迭代法求解是机械工程中求齿形系数和应力修正系数精度很高的一种方法，一次计算就可以求出一对相互齿合的齿轮的齿形系数和应力系数及所用迭代次数。方便实用

7 知识扩展

本节将在原有的不动点迭代法的基础之上进行推展到不动点迭代的加速的一些情况。

7.1 Steffensen 加速迭代

不动点迭代式（2.2）通常只有线性收敛，有时甚至不收敛，为加速迭代法的收敛性通常可以采用Steffensen 加速迭代。

设**()x x ?=是?的不动点，记*k k s x x =-，利用中值定理有：

**11**'

()()()k k k k k k k k s x x x x s x x x x

???ξξ++--===-- 在*x 与k x 之间通常'()k ?ξ依赖于k ，若'()x ?变化不大，设'()k C ?ξ≈于是有：

**1**21()()

k k k k x x C x x x x C x x +++-≈--≈-

从上式消去C ,则得

**21**

1k k k k x x x x x x x x

+++--≈--或***2

21()()()k k k x x x x x x ++--≈- 解得

2

2

*

212121()22k k k k k k k k k k k k x x x x x x x x x x x x x ++++++--≈=-

-+-+ 若记

2

121(),0,1,2......2k k k k k k

x x x x k x x x +++-=-=-+ （7.1）

用序列{}

k x 作为不动点*x 的新近似，一般说，它比迭代法（2.2）收敛更快，实际上迭代法式（2.1）（2.2）联立可改为：

2

1(),()

(),0,1,2......

2k k k k k k k k k k k

y x z y y x x x k z y x ??+==-=-=-+ （7.2） 称为Steffensen 迭代法，它是将来原不动点迭代（2.2）计算两次合并成一步得到，可改为另一种不动点迭代法：

1(),0,1,2......k k x x k ψ+== （7.3）

其中

2

[()]()(())2()k x x x x x x x

?ψ???-=-

-+ (7.4) 并有以下局部收敛定理：

若*x 为（7.4）定义的函数ψ的不动点，则*x 为?的不动点。反之，若*x 是?的不动点，设''*()x ?连续，且'*()1x ?≠，则*x 是ψ的不动点且迭代法（7.3）是二阶收敛的。

7.2 Aitken （埃特金）加速迭代法

设**()x x ?=是?的不动点，利用中值定理有：

**'*101010*

*

'

*

212121()()()()()()()()()()

x x x x x x x x x x x x x x x x ????ξ????ξ=?-=-=-=?-=-=-

若'()x ?变化不大，则可得：

''12()()?ξ?ξ≈

可推出：

*

*01**

21x x x x x x x x

--≈-- *1201210

()

2x x x x y x x x -≈-

=-+

即有如下图2所示关系关系

图2k x 与k y 的关系

进而可以得出其通项：

22

11221()()2k k k k k k k k k k

x x x y x x x x x x ++++-?=-=--+?

然后就有如下的收敛性：

*

1*

0lim k k

k y x x x +→∞-?-=超线性收敛

01234 . . .x x x x x 123

参考文献：

[1]John H.Mathews,Kurtis D.Fink.数值方法（第四版）[M].2010.4

[2]宜春学院学报405期.关于加速迭代的讨论.2009.3

[3] Gerald Recktenwald . 数值方法和MATLAB实现与应用[M].2004

[4] 金一庆 . 数值方法 [M].2011