浅谈牛顿迭代法在编程中的应用

牛顿迭代法(Newton's method)又称为牛顿-拉夫逊(拉弗森)方法

五次及以上多项式方程没有根式解(就是没有像二次方程那样的万能公式),这个是被伽罗瓦用群论做出的最著名的结论。

生产工作中还是有诸多求解高次方程或超越方程的真实需求(比如行星的轨道计算,往往就是涉及到很复杂的高次方程)

没有根式解不意味着方程解不出来,数学家也提供了很多方法,牛顿迭代法就是其中一种。

本文偏向应用

最经典的办法则是二分逼近,基于区间单调函数的0点存在定理

比如求 f(x)=\(\sqrt{x}=3\) 的解

先大概拟定一个根的范围[1,100]

通过f(50)=7.07 < 3 可以缩减区间为 [1,50]

通过f(25)=5 < 3 可以缩减区间为[1,25]

x=4 f(x)=2
x=10 f(x)=3.16

因为是单调函数则可以得到 x在区间[4,10]内

最终通过不断逼近得到答案

看完了二分法,我们再来看另一个快速求根的方法,和二分法一样,它也是迭代逼近的方法,但是逼近的速度更快。这个方法最早是牛顿提出的,因此也被称为牛顿迭代法,我想牛顿这个名字写出来,大家应该都能get到它的分量。

牛顿方法源于直觉,这种直觉本身有一定程度的合理性。

我们来看看收敛的充分条件:

若 f(x) 二阶可导,那么在待求的零点 周围存在一个区域,只要起始点 位于这个邻近区域内,那么牛顿-拉弗森方法必定收敛。

也就是说,在这个区域内,用切线代替曲线这个直觉是合理的。但是,因为我们不知道根点到底在哪里,所以起始点
选择就不一定在这个区域内,那么这个直觉就不可靠了。

1.随便确定一个数\(x_0\)
2.由导数定义可知 \(y-f(x_0)=f'(x_0)(x-x_0)\)
3.稍加计算 \(x_1=x_0-f(x_0)/f'(x_0)\)
既然是迭代,那么自然就有

double Newton(duoble x,int n){
     while(n--){
         x=x-f(x)/dx(x);
     }
     return x;
}

上面这个式子就是牛顿迭代法的迭代公式,这是一个非常牛的方法,比二分法要厉害得多,因为它的收敛速度更快,并且计算也并不复杂。

牛顿迭代法的名头看起来很唬人,但是原理真的不难,说白了只有一句话,就是通过切线去逼近,本质是借助泰勒级数,从初始值开始快速向零点逼近

不妨我们来尝试一下各种函数用牛顿迭代法的效果

1.最简单的 \(f(x)=x+1\)

\(f'(x)=1\)
\(x=x-(x+1)=-1\)
直接就解出 -1 ,只用了一次迭代
这种函数用于牛顿迭代法,确实大材小用

2.最简单的一种曲线 \(f(x)=x^2-4\)

\(f'(x)=2x\)

\(x=x-\frac{x^2-4}{2x}\)

\(x=x-\frac{x}{2}+\frac{2}{x}\)

\(x=\frac{x}{2}-\frac{2}{x}\)

就从100开始迭代,发现迭代了9次

50.020000
25.049984
12.604832
6.461085
3.540088
2.335002
2.024031
2.000143
2.000000

我们可以发现,逼近答案几乎是以一半一半的速度 ! !
我们再试一个数,这次大一些x=5000

2500.000400
1250.001000
625.002100
312.504250
156.258525
78.142062
39.096625
19.599468
9.901778
5.152873
2.964569
2.156919
2.005708
2.000008
2.000000

我们惊奇的发现,这个也是。
目前我们可以暂且认定二次函数复杂度是 \(O(log x)\) (x指的是与答案的差)
不过我们现在得到的都是其中一个解,我们把x=-1000,就可以去求一下负解
对于这种两解方程,可以设置一个正无穷,一个负无穷

3. \(f(x)=ln(x)-1\)

这个和之前的都不同
解这个方程必须要注意,因为定义域问题 \(x>0\)
如果初始x与答案差别过大,很容易出现报错
比如:\(x=1000\)

-4907.755279

爆了一个这个东西
因为对数运算的上标和下表必须>0
而之前的定义域都是R,不存在这样的问题
怎么办呢?其实很简单

double Newton(duoble x,int n){
     while(n--){
         x=max(0.01,x);
         x=x-f(x)/dx(x);
     }
     return x;
}

整个过程保证在定义域内就好了

  1. 一个超越方程 cosx+x=1
    这次我就省略求导过程了
    这个带三角函数,会怎么样呢?
    从x=1000开始尝试

-4773.796940
438522.335367
7235.249930
1006.875512
-82651695.047529
744885946.410124
-844300663.013884
160826680.451773
-1313204400.701244
7428745450.461179
-8159334547.905384
27158929610.750668
11670359502.680054
754141218.371896
364818504.965956
-557131947.826826
-129458645.156508
-62926439.507690
91972060.494912
-46017074.491629
118586484.500184
-1693025646.669981
-406239679.926741
-177004647.243325
100583519.521505
-812016349.392783
1866388963.645748
-42400283650.161812
-13600371300.391312
-6193017186.739635
5425058253.690373
2334777092.105400
-208560177707.866058
-97224777319.680511
-27363005251.129093
18767623801.335831
-14475614337.709970
-6270820388.583565
-2108067808.927387
-972584844.203155
235203161.752808
111045994.146395
45312266.220592
-199964166.191071
-11192194.933260
-4975516.590042
-1334745.837151
-500082.021570
95737.701323
-248735.015858
-88270.537451
1272760.735764
-7358582.063434
221629465.505894
34899989.085697
17325319.061259
-120185094.177484
1002640851.882419
497684164.620903
186246003.982280
-1147881086.324440
262233347.990792
82756105.034660
-997521893.542680
-486098632.306040
-223350507.882394
-93660535.654878
-46468861.942258
-14868377.712947
205688129.226440
2189543.868569
1090351.500038
508349.917610
-692014.197815
-76142.614860
-10811.602831
600993.017739
-38283.545295
-3387.208519
-1189.062481
-593.909468
105.429084
52.668515
-103.184560
-31.654279
-6.023266
2.128651
-1.824217
-0.261799
-0.026759
-0.000349
-0.000000

一共92次迭代,在三角函数内部的复杂度并不是 \(O(log x)\)
根据我个人经验:带三角函数的方程,速率都不太稳定
而且时常带有多解的情况
所以请尽量缩短值域

无法收敛的情况

令人遗憾的是并不是所有方程使用牛顿迭代法都可以有这么好的效果,对于一些方程,甚至可能会出现越走越偏的情况。我们再举个例子,比如方程f(x)=\(x^\frac{1}{3}\)。如果我们画出它的迭代过程,是这样的:

此处根点很显然是0点,而 \(f'(x)\) 是不存在的。

具体多近才能保证收敛?一般来说是无法确定的,应用时往往会尽量选取一个非常靠近的值。

牛顿迭代法最常用的还是计算机求解开平方根,可以容易找到近似值

最后总结

应用牛顿方法,要注意以下问题:

函数在整个定义域内最好是二阶可导的

起始点对求根计算影响重大,可以增加一些别的判断手段进行试错

本文部分内容参考:

《用计算机求二类收敛级数和的近似值》许致一

posted @ 2020-04-17 15:10  白木偶君  阅读(383)  评论(0编辑  收藏  举报