数值方法笔记2：解决非线性方程

1. 不动点定理及其条件验证

不动点定义：$P=g(P)$
不定点迭代：$x_{k+1}=g\left(x_k\right)$
定理： 如果$g(x_k)$是连续的并且序列$x_k$是收敛的，$x_k$收敛到方程的解：$x=g(x)$
$$x^{\ast }=g\left(x^{\ast }\right)\text{ and }{x}_k->x^{\ast }$$

定理： 假设
(1) 对于$g(x)$，$g^{\prime }(x)\in C[a,b]$（连续）
(2) $K$是一个正的常数
(3) $p_0\in (a,b)$
(4) $g(x)\in [a,b],\forall x\in [a,b]$
那么
(a) 如果$\left|{\mathrm{g}}^{\prime }(x)\right|\le \mathrm{K}<1，\forall x\in [\mathrm{a},\mathrm{b}],{\mathrm{x}}_{\mathrm{k}+1}=\mathrm{g}\left({\mathrm{x}}_{\mathrm{k}}\right)$收敛。
(b) 如果$\left|{\mathrm{g}}^{\prime }(x)\right|>1，\forall x\in [\mathrm{a},\mathrm{b}],{\mathrm{x}}_{\mathrm{k}+1}=\mathrm{g}\left({\mathrm{x}}_{\mathrm{k}}\right)$不收敛

曲线的切线斜率$k\in (-1,1)$看下面的图逐渐收敛：

曲线的切线斜率$k\in [-\infty ,-1)\cup (1,\infty ]$看下面的图不收敛：

综上所述，不动点迭代满足最重要的是：

$\begin{array}{rlr}& (1)\left|{\mathrm{g}}^{\prime }(x)\right|\le \mathrm{K}<1,\forall x\in [\mathrm{a},\mathrm{b}]& 【g^{\prime }(x)的边界条件】\\ & (2)\mathrm{g}(x)\in [\mathrm{a},\mathrm{b}]，\forall x\in [\mathrm{a},\mathrm{b}]，并且有g([a,b])\subset [a,b]& 【g(x)的边界条件】\end{array}$

单调和非单调要分别判断边界条件，单调的g(x)的范围看端点就可以了，非单调还要看极值点。

2. 收敛阶、收敛检测与收敛加速

定义;
$$\left|x_{k+1}-x^{\ast }\right|\le C{\left|x_k-x^{\ast }\right|}^p,k>M\text{, for }C>0,p>0$$
或
$$\underset{k\to \infty }{\lim }\frac{\left|x_{k+1}-x^{\ast }\right|}{{\left|x_k-x^{\ast }\right|}^p}=C$$
为$p$阶收敛
其中：
$p=1$ , 线性收敛（linear convergence）
$1<p<2$ , 超线性收敛（superlinear convergence）
$p=2$ , 平方收敛（square convergence）

2.1 如何估计不动点迭代的收敛阶$x_{k+1}=g\left(x_k\right)$

定理;设$x^{\ast }$是最优解，如果$g^{\prime }\left(x^{\ast }\right)=g^{\prime \prime }\left(x^{\ast }\right)=\dots =g^{(p-1)}\left(x^{\ast }\right)=0$, $g^{(p)}\left(x^{\ast }\right)\ne 0$， $x_{k+1}=g\left(x_k\right)$是 $p$阶收敛。

证明：
$$\begin{array}{rl}{x}_{k+1}& =g\left(x_k\right)=g\left(x^{\ast }\right)+g^{\prime }\left(x^{\ast }\right)\left(x_k-x^{\ast }\right)+\dots \\ & +\frac{g^{(p-1)}\left(x^{\ast }\right){\left(x_k-x^{\ast }\right)}^{p-1}}{(p-1)!}+\frac{g^{(p)}(\xi ){\left(x_k-x^{\ast }\right)}^p}{p!},【\xi \in \left[x_k,x^{\ast }\right]或\left[x^{\ast },x_k\right]】\\ \Rightarrow & x_{k+1}=x^{\ast }+\frac{g^{(p)}(\xi ){\left(x_k-x^{\ast }\right)}^p}{p!}\\ \Rightarrow & \frac{x_{k+1}-x^{\ast }}{{\left(x_k-x^{\ast }\right)}^p}=\frac{g^{(p)}(\xi )}{p!}\to \frac{g^{(p)}\left(x^{\ast }\right)}{p!}\end{array}$$

2.2 给定精度的情况下，如何预测不动点迭代需要迭代的次数

定义 L = max ⁡ x ∈ [ a , b ] { ∣ g ′ ( x ) ∣ } < 1 L=\max _{x \in[a, b]}\left\{\left|g^{\prime}(x)\right|\right\}<1 L=maxx∈[a,b]​{∣g′(x)∣}<1
迭代的次数满足：$k\ge \ln \left(\epsilon (1-L)/\left|x_1-x_0\right|\right)/\ln L$

证明：
$$\begin{array}{r}{x}_{k+1}=g\left(x_k\right)=g\left(x^{\ast }\right)+g^{\prime }(\xi )\left(x_k-x^{\ast }\right)=x^{\ast }+g^{\prime }(\xi )\left(x_k-x^{\ast }\right)\\ \Rightarrow \left|x_{k+1}-x^{\ast }\right|\le \left|g^{\prime }(\xi )\Vert \left(x_k-x^{\ast }\right)\right|\le L\left|x_k-x^{\ast }\right|\le L^k\left|x_1-x^{\ast }\right|\end{array}$$
又有
$$\left|x_{k+1}-x_k\right|=\left|g\left(x_k\right)-g\left(x_{k-1}\right)\right|\le L\left|x_k-x_{k-1}\right|\le L^k\left|x_1-x_0\right|$$
于是有：
$$\begin{array}{rl}& \left|x_{k+q}-x_k\right|\le \left|x_{k+q}-x_{k+q-1}\right|+\left|x_{k+q-1}-x_{k+q-2}\right|+\dots +\left|x_{k+1}-x_k\right|\\ & \le \left(L^{q-1}+L^{q-2}+\dots +1\right)\left|x_{k+1}-x_k\right|\\ & <\left(1+L+L^2+\dots +L^{q-1}+\dots \right)\left|x_{k+1}-x_k\right|\\ & =\frac{1}{1-L}\left|x_{k+1}-x_k\right|\\ & \le \frac{L^k}{1-L}\left|x_1-x_0\right|\end{array}$$

让$q\to \infty$有

$$\left|x^{\ast }-x_k\right|\le \frac{1}{1-L}\left|x_{k+1}-x_k\right|\le \frac{L^k}{1-L}\left|x_1-x_0\right|$$

于是：

$$\frac{L^k}{1-L}\left|x_1-x_0\right|\le \epsilon \Rightarrow k\ge \ln \left(\epsilon (1-L)/\left|x_1-x_0\right|\right)/\ln L$$

2.3 如何加快收敛的速度

$$\begin{array}{rl}& x_{k+1}-x^{\ast }\approx L\left(x_k-x^{\ast }\right)\\ & x_{k+2}-x^{\ast }\approx L\left(x_{k+1}-x^{\ast }\right)\\ & \frac{x_{k+1}-x^{\ast }}{x_{k+2}-x^{\ast }}\approx \frac{x_k-x^{\ast }}{x_{k+1}-x^{\ast }}\Rightarrow x^{\ast }\approx x_k-\frac{{\left(x_{k+1}-x_k\right)}^2}{x_{k+2}-2x_{k+1}+x_k}=x^{\Delta }\end{array}$$

根据上面的思路我们可以：

$$\begin{array}{rl}Iteration& {\bar{x}}_{k+1}=g\left(x_k\right)\\ Onemore& {\hat{x}}_{k+1}=g\left({\bar{x}}_{k+1}\right)\\ Tospeedup& x_{k+1}=x_k-\frac{{\left({\bar{x}}_{k+1}-x_k\right)}^2}{{\hat{x}}_{k+1}-2{\bar{x}}_{k+1}+x_k}\end{array}$$

2.4 停止不定点迭代的条件

当 L = max ⁡ x ∈ [ a , b ] { ∣ g ′ ( x ) ∣ } < 1 L=\max _{x \in[a, b]}\left\{\left|g^{\prime}(x)\right|\right\}<1 L=maxx∈[a,b]​{∣g′(x)∣}<1时，可以使用下面的条件：

$$\left|x_{\mathrm{k}+1}-x_{\mathrm{k}}\right|<\mathrm{eps}$$

2.5 不动点迭代的两个缺点

1. 很难估计 L ( max ⁡ x ∈ [ a , b ] { ∣ g ′ ( x ) ∣ } ) L(\max _{x \in[a, b]}\left\{\left|g^{\prime}(x)\right|\right\}) L(maxx∈[a,b]​{∣g′(x)∣})
2. $L<1$时无法收敛。

3. 应用：如何求解非线性方程组$f(x)=0$的解

3.1 二分法(Bisection Method of Bolzano)

算法的流程：

1. 用一个区间找到一个根。
2. 用中点分割该区间。
3. 选择其中的一个子区间作为新的位置。

$$\begin{array}{rl}& a=x_0,b=x_0+h\\ & c=\frac{a+b}{2}\\ & f(a)f(b)<0,\end{array}$$
于是：

$$\begin{array}{rl}& [a,b]\to \left[a_1,b_1\right]\to \left[a_2,b_2\right]\to \dots \to \left[a_n,b_n\right]\\ & a=a_0\le a_1\le \cdots \le a_n\le \cdots \le r\le \cdots \le b_n\le \cdots \le b_1\le b_0=b\end{array}$$
定义$r$是精确解。
$$\begin{array}{rl}& \left|r-c_n\right|\le \frac{b-a}{2^{n+1}},\text{ for }n=0,1,2,\dots \\ & c_n=\frac{a_n+b_n}{2}\end{array}$$

迭代次数N：

$$\begin{array}{rl}& \left|r-c_n\right|\le \frac{b-a}{2^{n+1}}<\delta \\ & 2^{n+1}>\frac{b-a}{\delta }\\ & (n+1)\ln 2>\ln (b-a)-\ln \delta \\ & n+1>\frac{\ln (b-a)-\ln \delta }{\ln 2}\\ & N=\mathrm{int}\left(\frac{\ln (b-a)-\ln \delta }{\ln 2}\right)\end{array}$$

简单地利用二分法可以判断区间内有没有零点（区间内有变号【可取最大值和最小值】）

3.2 试位法(False Position Method)

算法的流程：

1. 用一个区间找到一个根。
2. 以割线与X轴的交点划分区间。(过程中仍然保证端点的异号，让区间包含零点)
3. 选择其中一个子区间作为新的位置。

$$\begin{gathered} c=b-\frac{f(b)(b-a)}{f(b)-f(a)} \\ {c}_1\to c_2\to \dots \to r \\ \left[a_n,b_n\right]\to [a,c]:=\left[a_{n+1},b_{n+1}\right] \end{gathered}$$

缺点：在凹函数下不适用，不会收敛。

3.3 牛顿-拉夫逊方法(Newton-Raphson method)

我们知道不动点迭代，能不能用到求解非线性方程组呢？

使用泰勒展式：

$$f(x_{k+1})=f\left(x_k\right)+f^{\prime }\left(x_k\right)(x_{k+1}-x_k)+O\left(|d{|}^2\right)=0$$

于是我们可以让

$$f\left(x_{\mathrm{k}}\right)+f^{\prime }\left(x_{\mathrm{k}}\right)\left(x_{k+1}-x_k\right)=0$$

使得：

$$x_{\mathrm{k}+1}=x_k-f\left(x_{\mathrm{k}}\right)/f^{\prime }\left(x_{\mathrm{k}}\right)=g(x_k)$$

总结Newton-Raphson方法即：

$$\begin{array}{l}f(x)=0\\ x=g(x)=x-\frac{f(x)}{f^{\prime }(x)}\\ x_{k+1}=g\left(x_k\right)=x_k-\frac{f\left(x_k\right)}{f^{\prime }\left(x_k\right)}\end{array}$$

我们可以证明在解的附近，Newton-Raphson方法是收敛的。

证明：

$$g(x)=x-f(x)/f^{\prime }(x)$$

$$g^{\prime }(x)=1-\frac{f^{\prime }(x)f^{\prime }(x)-f(x)f^{\prime \prime }(x)}{{\left[f^{\prime }(x)\right]}^2}=\frac{f(x)f^{\prime \prime }(x)}{{\left[f^{\prime }(x)\right]}^2}$$

我们知道不动点的条件是$\left|g^{\prime }(x)\right|<K<1$，当我们取的邻域足够小，条件$g([a,b])\subset [a,b]$会满足，注意到$f(x^{\ast })=0$，在解的邻域附近，因为$f(x)=0$，所以$g^{\prime }(x)=0$。

各种条件下的推导（不做要求，想了解可以看一下）

1. $f^{\prime }\left(x^{\ast }\right)>0\text{ and }{f}^{\prime \prime }\left(x^{\ast }\right)<0,g\left(\left[x^{\ast }-\delta ,x^{\ast }+\delta \right]\right)\subset \left[x^{\ast }-\delta ,x^{\ast }+\delta \right]$
   
   $$\begin{array}{rl}& x^{\ast }-\delta <g\left(x^{\ast }-\delta \right)=\left(x^{\ast }-\delta \right)-\frac{f\left(x^{\ast }-\delta \right)}{f^{\prime }\left(x^{\ast }-\delta \right)}\\ \iff & 0<-\frac{f\left(x^{\ast }-\delta \right)}{f^{\prime }\left(x^{\ast }-\delta \right)}\\ \iff & \frac{f\left(x^{\ast }-\delta \right)}{f^{\prime }\left(x^{\ast }-\delta \right)}<0\\ \iff & f\left(x^{\ast }-\delta \right)<0\\ \iff & f\left(x^{\ast }\right)-f^{\prime }(\xi )\delta <0【\xi \in [x^{\ast }-\delta ,x^{\ast }]】\\ \iff & -f^{\prime }(\xi )\delta <0\\ \Rightarrow & \exists {\delta }_1>0,f^{\prime }(\xi )>0,\text{ for }{x}^{\ast }-\xi <{\delta }_1\end{array}$$
   又有
   $$\begin{array}{rl}& f^{\prime \prime }\left(x^{\ast }\right)<0\\ \Rightarrow & \exists {\delta }_2>0,f^{\prime \prime }(x)<0【保号性】\\ \Rightarrow & g^{\prime }(x)=\frac{f(x)f^{\prime \prime }(x)}{{\left[f^{\prime }(x)\right]}^2}>0,\text{ for }{x}^{\ast }-x<{\delta }_2\\ & 【{\delta }_2足够小，导数保号性，f^{\prime }(x)>0,x<x^{\ast },f(x^{\ast })=0，f(x)<0】\end{array}$$
   当 δ < min ⁡ { δ 1 , δ 2 } \delta<\min\{\delta_1,\delta_2\} δ<min{δ1​,δ2​}有：
   $$x^{\ast }-\delta <g\left(x^{\ast }-\delta \right)<g(x),\text{ for }{x}^{\ast }-x<\delta$$
2. $f^{\prime }\left(x^{\ast }\right)>0\text{ and }{f}^{\prime \prime }\left(x^{\ast }\right)<0,g\left(\left[x^{\ast }-\delta ,x^{\ast }+\delta \right]\right)\subset \left[x^{\ast }-\delta ,x^{\ast }+\delta \right]$
   
   $$\begin{array}{rl}& x^{\ast }-\delta <g\left(x^{\ast }-\delta \right)=\left(x^{\ast }-\delta \right)-\frac{f\left(x^{\ast }-\delta \right)}{f^{\prime }\left(x^{\ast }-\delta \right)}\\ \iff & 0<-\frac{f\left(x^{\ast }-\delta \right)}{f^{\prime }\left(x^{\ast }-\delta \right)}\\ \iff & \frac{f\left(x^{\ast }-\delta \right)}{f^{\prime }\left(x^{\ast }-\delta \right)}<0\\ \iff & f\left(x^{\ast }-\delta \right)<0\\ \iff & f\left(x^{\ast }\right)-f^{\prime }(\xi )\delta <0【\xi \in [x^{\ast }-\delta ,x^{\ast }]】\\ \iff & -f^{\prime }(\xi )\delta <0\\ \Rightarrow & \exists {\delta }_1>0,f^{\prime }(\xi )>0,\text{ for }{x}^{\ast }-\xi <{\delta }_1\end{array}$$
   又有
   $$\begin{array}{rl}& f^{\prime \prime }\left(x^{\ast }\right)>0\\ \Rightarrow & \exists {\delta }_2>0,f^{\prime \prime }(x)<0【保号性】\\ \Rightarrow & g^{\prime }(x)=\frac{f(x)f^{\prime \prime }(x)}{{\left[f^{\prime }(x)\right]}^2}<0,\text{ for }{x}^{\ast }-x<{\delta }_2\\ & 【{\delta }_2足够小，导数保号性，f^{\prime }(x)>0,x<x^{\ast },f(x^{\ast })=0，f(x)<0】\end{array}$$
   当 δ < min ⁡ { δ 1 , δ 2 } \delta<\min\{\delta_1,\delta_2\} δ<min{δ1​,δ2​}有：
   $$x^{\ast }-\delta <x^{\ast }=g\left(x^{\ast }\right)<g(x),\text{ for }{x}^{\ast }-x<\delta ,x<x^{\ast }$$

注意Newton-Raphson方法对于单根是二阶收敛(二次收敛)【quadratic convergence】

$$\left|E_{n+1}\right|\approx \frac{\left|f^{\prime \prime }(p)\right|}{2\left|f^{\prime }(p)\right|}{\left|E_n\right|}^{2}n\to \infty$$

证明：

而对于多重根是线性（一次）收敛，收敛速度降低。

$$\left|E_{n+1}\right|\approx \frac{M-1}{M}\left|E_n\right|n\to \infty$$

证明：

如果出现了多重根$p^{\ast }$，我们看到在$f^{\prime }(p^{\ast })=0$，Newton-Raphson方法的分母会出现0.然而一般来说，分子$f(p_k)$要比分母$f^{\prime }(p_k)$先出现0，所以Newton-Raphson方法一般还是可以用的。

Newton-Raphson方法的问题：

1.分母可能为0，除以零是不允许的。
2.收敛到一个不同的根，或发散。
3.产生一个循环序列。
4.产生一个发散的振荡序列。

由于多重根线性收敛的问题，可以考虑Newton-Raphson方法加速：
$$p_k=p_{k-1}-\frac{Mf\left(p_{k-1}\right)}{f^{\prime }\left(p_{k-1}\right)}M>1$$

证明：

3.4 割线法(Secant Method)

当Newton-Raphson的导数不好显式表达的时候，可以通过两端点的直线的斜率来近似导数。

我们有：

$$x_{k+2}=g\left(x_k,x_{k+1}\right)=x_{k+1}-\frac{f\left(x_{k+1}\right)\left(x_{k+1}-x_k\right)}{f\left(x_{k+1}\right)-f\left(x_k\right)}$$

3.5 Aitken过程加速

使用不定点的迭代，Aitken过程加速又称为史蒂芬森加速（Steffensen’s acceleration）.注意，只对一阶方法有效。

$$\underset{n\to \infty }{\lim }\frac{p-p_{n+1}}{p-p_n}=A,p\approx \frac{p_{n+2}{p}_n-p_{n+1}^2}{p_{n+2}-2p_{n+1}+p_n}=q_n$$

3.6 Muller方法（Muller’s method）

给定三个初始值$\left(p_0,f\left(p_0\right)\right),\left(p_1,f\left(p_1\right)\right),\left(p_2,f\left(p_2\right)\right)$

令
$$\begin{array}{rl}& t=x-p_2\\ & h_0=p_0-p_2,h_1=p_1-p_2\end{array}$$

我们使用二次函数计算下一个点：

$$y=a{t}^2+bt+c$$

则有：

$$\begin{array}{rl}t=h_0:a{h}_0^2+b{h}_0+c=f_0& \Rightarrow a{h}_0^2+b{h}_0=f_0-c=e_0\\ t=h_1:a{h}_1^2+b{h}_1+c=f_1& \Rightarrow a{h}_1^2+b{h}_1=f_1-c=e_1\\ t=0:a{0}^2+b0+c=f_2& \Rightarrow c=f_2\end{array}$$

解得：

$$a=\frac{e_0{h}_1-e_1{h}_0}{h_1{h}_0^2-h_0{h}_1^2},b=\frac{e_1{h}_0^2-e_0{h}_1^2}{h_1{h}_0^2-h_0{h}_1^2}$$

于是得到：

a t 2 + b t + c = 0 : t = z 1 , z 2 ⇒ z i = − 2 c b ± b 2 − 4 a c z = arg ⁡ min ⁡ { ∣ z i ∣ } 【对于一个复数，在计算中只保留其实数部分】 \begin{aligned} &a t^{2}+b t+c=0: \quad t=z_{1}, z_{2} \Rightarrow z_{i}=\frac{-2 c}{b \pm \sqrt{b^{2}-4 a c}} \\ &z=\arg \min \left\{\left|z_{i}\right|\right\}【\text{对于一个复数，在计算中只保留其实数部分}】 \end{aligned} ​at2+bt+c=0:t=z1​,z2​⇒zi​=b±b2−4ac ​−2c​z=argmin{∣zi​∣}【对于一个复数，在计算中只保留其实数部分】​

$$p_3=p_2+z$$

继续得到$\left({\bar{p}}_1,{\bar{p}}_2,p_3\right)$，其中${\bar{p}}_1,{\bar{p}}_2$是距离$p_3$最近的两个点。

4. 其他问题

4.1 如何寻找初值

例如

可以有两个判断条件：

1. 【针对$r_1$和$r_2$】
   $$f\left(x_{k-1}\right)f\left(x_k\right)<0[a,b]=\left[x_{k-1},x_k\right]$$

2. 【针对$r_3$】
   $$\left|f\left(x_k\right)\right|<\epsilon \text{ 并且}\left(f\left(x_k\right)-f\left(x_{k-1}\right)\right)\left(f\left(x_{k+1}\right)-f\left(x_k\right)\right)<0[a,b]=\left[x_{k-1},x_{k+1}\right]$$

4.2 收敛条件

可以有两个收敛条件：

1. 根据纵坐标
$$\left|f\left(x_k\right)\right|<\epsilon$$

误差为：${\text{Error}}_x=\left|x_k-r\right|$

2. 根据横坐标

$$\left|x_k-x_{k-1}\right|<\delta$$

由以下推出：

$$\left|x_k-r\right|<\delta \Rightarrow \left|x_k-x_{k-1}\right|<\delta$$

误差为：  Error  f = max ⁡ { ∣ f ( r − δ ) ∣ , ∣ f ( r + δ ) ∣ } \text { Error }_{f}=\max \{|f(r-\delta)|,|f(r+\delta)|\}  Error f​=max{∣f(r−δ)∣,∣f(r+δ)∣}

3. 我们也可以把上面两个进行组合：

$$\left|f\left(x_k\right)\right|<\epsilon \text{并且}\left|x_k-r\right|<\delta$$

4. 如果针对Newton-Raphson问题，我们还可以有如下的判断标准：

①$f^{\prime }(r)\ne 0$
②$x_0\in [r-\delta ,r+\delta ]$, $\delta$足够小。

4.3 算法的收敛速度对比

4.4 算法的选择

单根：
Newton-Raphson方法

双根（当分母为0失效）：
Newton-Raphson方法
Steffensen’s method