【转】布伦特方法(Brent‘s method)---求根方法

1. 解方程(Solving Equations)

1.1 二分法(The Bisection Method)

定义1.1 对于方程

f(x)

,如果有

f(r) = 0

,则说 

x = r

 是

f(x)

 的一个根。

定理1.2 令 

f

 是 

[a, b]

 上的连续函数,满足 

f(a)f(b) < 0

 , 则 

f

 在 

a

 和 

b

 之间存在一个根。

二分法算法流程:

二分法误差:

在区间 

[a, b]

n

 次二分之后,得到最终区间 

[a_{n}, b_{n}]

 其长度为 

(b - a) / 2^{n}

 。选取中间点 

x_{c} = (a_{n} + b_{n}) / 2

最为根 

r

 的最有估计,则误差为:

定义1.3 如果误差小于 

0.5 \times 10^{-p}

,则称其精确到小数点后 

p

 位。

二分法迭代次数:

只要确定了期望的精度,就可以依据二分法误差计算公式得到所需要的迭代次数。

1.2 不动点迭代(Fixed-Point Iteration)

定义1.4 对于方程 

g

,存在 

r

 ,使得 

g(r) = r

,就称 

r

 是方程 

g

 的不动点。

例如:

FPI 算法流程:

当方程可写作 

g(x) = x

 的形式时,则使用初值 

x_{0}

 进行 FPI:

FPI 可能收敛,也可能不收敛,如果收敛到 

r

,则 

r

 就是 

g

 的不动点,也是 

f(x) = g(x) - x = 0

 的解

不同的 

x = g(x)

 形式

同一个方程写成不同的 

x = g(x)

,由于 

g(x)

 的不同,其能不能收敛,收敛快慢也不同:

:对于 

x^3 + x - 1 = 0

,其 

g(x)

 可写成:

g(x) = 1 - x^3

                 (a)

g(x) = \sqrt[3]{1 - x}

               (b) 

g(x) = \frac{1 + 2x^3}{1 + 3x^2}

               (c)

令 

x_{0} = 0.5

,分别进行 FPI,其每次迭代数值,及 cobweb 图 :

(a)  不收敛

(b) 收敛较慢

       

(c) 快速收敛

      

不动点迭代(FPI)的线性收敛

通过对 

g(x)

 几何观察,其不动点 

r

 附近的斜率的绝对值,即 

|{g}'(\tilde{r})|

 如果 小于 1,则FPI收敛,大于 1,则FPI不收敛

定义1.5 令 

e_{i}

 表示FPI过程中第 

i

 步的误差,如果:

满足线性收敛,收敛率为

S

注: 关于 线性收敛,超线性收敛,r阶收敛以及 Rate of convergence令见下文 定义1.10。

定理1.6 设 

g

 连续可微,

g(r) = r

, 且 

S = |g'(r)| < 1

 。那么就称对于一个足够接近 

r

 的初值估计,FPI以速度 

S

 线性收敛到不动点

r

.

定义1.7 如果一个迭代方法,能够以一个足够接近 

r

 的初值估计收敛到 

r

,则称该方法局部收敛(Locally Convergent)

:使用FPI计算 

\sqrt{2}

x = \sqrt{2} \Rightarrow x^2 - 2 = 0 \Rightarrow x + \frac{2}{x} -2x =0 \Rightarrow \frac{x + \frac{2}{x}}{2} = x

,即求 

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

 的不动点,运用FPI,收敛到 1.414215686

分析此时 

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

 导数 

g' = \frac{1}{2} - x^{-2}

,收敛率

S = |\frac{1}{2} - (\sqrt{2}^{-2})| = 0

,收敛,且收敛很快。

终止条件(Stopping Criteria)

不同于二分法,FPI无固定的误差公式,无法预测收敛所需步数,故而需要给出终止条件

(1.18)中 

\theta > 0

 通常用于解在 0 附近的情况。

另外,好的FPI(代码)实现通常要设置最大迭代次数,以防止收敛失败

二分法和FPI比较

二分法FPI
收敛性保证线性收敛局部收敛时,则线性收敛
每步计算每步只需一次函数求值一次
速度每步去掉1/2不确定性不确定性会乘上

S=|g'(r)|

,故而可能更慢或更快,要看

S

比1/2大还是小
其他牛顿方法是FPI的改善,其 

S = 0

1.3 精度(Accuracy)

前向误差和后向误差 (Forward error & Backward error)

定义1.8 设 

r

是方程 

f

 的一个根,即满足 

f(r) = 0

, 

x_{a}

 是根 

r

 的近似值,对于求解根的问题,该近似解 

x_{a}

 的前向误差是 

|r - x_{a}|

后向误差是 

|f(x_{a})|

 。

关于“前向”和“后向”:

Backward error is on the left or input (problem data) side. It is the amount we would need to change the problem (the function f ) to make the equation balance with the output approximation 

x_{a}

. This amount is 

|f(x_{a})|

Forward error is the error on the right or output (problem solution) side. It is the amount we would need to change the approximate solution
to make it correct, which is 

|r - x_{a}|

.”

思考:对于方程 

f(x) = 0

 有根 

r

,写成 

g(x) = x

 的形式,即有 

g(r) = r

,设

x_{a}

 为近似值,则对应的:

前向误差

|r - x_{a}|

后向误差

|x_{a} - g(x_{a})|

相对前向误差(relative forward error)

\frac{|r - x_{a}|}{|r|}

相对后向误差(relative backward error)

\frac{|x_{a} - g(x_{a})|}{|x_{a}|}

定义1.9 设 

r

 是可微函数 

f

 的 一个根,即满足 

f(r) = 0

,如果有:

 

就称 

f

 在 

r

 处有 

m

 重数(multiplicity m)。如果重数大于1,则称 

f

 在 

r

 处有重根(multiple root),当重数为1时,称为单根(simple)。

Wilkinson多项式

由式(1.19)可知根为1 到 20 的实数。但是若是给出的是式(1.20),使用迭代方法求解,“its evaluation suffers from cancellation of nearly equal, large numbers.”  小误差导致最终的大误差。

根求解的敏感性(Sensitivity of root-finding)

当输入小的误差,导致求解后大的输出误差时,这一情况被称作敏感性

根的敏感性公式:

误差放大因子(error magnification factor):

即: 相对前向误差 / 相对后向误差

条件数(condition number)

"The condition number of a problem is defined to be the maximum error magnification over all input changes, or at least all changes of a prescribed type."  条件数是误差放大的一总度量。

下图来自 维基百科-条件数

思考条件数 也可以看作是 最大误差放大因子

条件数为1左右的问题称为 良态问题(well-conditioned),条件数高的问题称为 病态问题(ill-conditioned)。

1.4 牛顿法(Newton's Method)

牛顿-拉弗逊法(Newton-Raphson Method),通常比之前提的线性收敛的收敛速度更快。

牛顿法算法流程:

牛顿法的二次收敛(Quadratic convergence of Newton’s Method)

定义1.10 令 

e_{i}

 为某迭代法第 

i

 步后的误差,则满足:

,称该方法二次收敛

定理1.11 设 

f

 是二阶连续可微函数,且 

f(r) = 0

 ,若 

f'(r) \neq 0

,则称牛顿法局部二次收敛到 

r

。第 

i

 步的误差为满足:

证明:

牛顿法的线性收敛(Linear convergence of Newton’s Method)

定理1.11 并不是说牛顿法总是能够二次收敛。

定理1.12 设在区间 

[a, b]

 上,有 

(m + 1)

 阶连续可微函数 

f

,其 

r

 有 

m

 重根,那么牛顿法局部收敛到 

r

,且 

i

 步的误差 

e_{i}

 满足:

,其中 

S = (m - 1) / m

注:依据定理1.11 和定理1.12,可以知道牛顿法在 

r

 是单根情况下为二次收敛

r

 为重根时,牛顿法的收敛率同二分法和FPI一样同处于线性收敛

改进牛顿法(Modified Newton’s Method)

定理1.13 如果 

f

 是在区间 

[a, b]

 上的 

(m + 1)

 阶连续可微函数,其有 

m

 重根 

r

,那么牛顿法改进为: 

局部二次收敛到 

r

 。

1.5 无导数的根求解

除了多重根的情况,牛顿法之所以比二分法和FPI收敛更快,是因为其运用了更多的信息,就是导数。但有的时候导数难以计算,所以有以下无需导数的迭代求解方法:

割线法(Secant Method)

同牛顿法相比,即用差商替代导数。且需要两个初始估值。并以超线性(superlinear)收敛到 

r

,收敛速度介于线性和二次收敛之间。割线法有三个变体

变体一·试位法(False Position 或 Regula Falsi)

与二分法有些类似,中点由类割线方法替代

注:试位法初始表现出对二分法和割线法的性能提升,因其融合了二者的优点,但二分法能够保证每次去掉1/2的不确定性,而试位法并不能有此保证,某些时候可能收敛的非常慢:

变体二·穆勒法(Muller's Method)

使用前三个迭代值

x_{0}, x_{1}, x_{2}

画出一个抛物线 

y = p(x)

,该抛物线与x轴相交,通常会有0到2个交点。如果有两个交点,那么距离

x_{2}

的那个交点作为

x_{3}

。如果不相交,则会有复数解。

变体三·逆二次插值(Inverse Quadratic Interpolation 即IQI)

使用形如 

x = p(y)

 的抛物线替代穆勒法中的 

y = p(x)

 。那么该抛物线与x轴只有一个交点。

经过 

(a, A), (b, B), (c, C)

 三个点的二次多项式 

y = p(x)

 为:

这是拉格朗日插值(Lagrange interpolation)的一个例子(后头讲插值时会细讲)。

穆勒法和IQI比较示意:

穆勒法和IQI均比割线法收敛快,这是因为使用的更高阶的插值。

Brent法

一种混合方法,其融合了之前出现的一些方法。试想结合了二分法一定收敛,且其他复杂方法的快速收敛特性,岂不是能爽到!

同样,Brent法运用于连续函数 

f

 ,边界是

a, b

,且

f(a)(b) < 0

Brent法记录当前点 

x_{i}

,其有最优后向误差

x_{i}

 that is best in the sense of backward error),且有包含有根的区间 

[a_{i}, b_{i}]

 。总督来说就是:如果满足1)后向误差得到改进(improve)且 2)包含根的区间

[a_{i}, b_{i}]

至少减半,那么尝试IQI法取更新 

x_{i}, a_{i}, b_{i}

 的值;如果不满足以上两点,则使用割线法更新

x_{i}, a_{i}, b_{i}

;如果还是失败了,就用二分法,这样能保证每步迭代以至少去除一半不确定性。

版权声明:本文为博主原创文章,遵循 CC 4.0 BY-SA 版权协议,转载请附上原文出处链接和本声明。

本文链接:数值分析一下咯(一)_shanexn的博客-CSDN博客_wilkinson多项式

  • 2
    点赞
  • 15
    收藏
    觉得还不错? 一键收藏
  • 1
    评论
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值