非线性方程的解
目录From 2022-12-2 To 2022-12-
Learning from 物理学中的非线性方程的逐步搜索法和二分法 求解非线性方程的迭代法 弦截法
带有空气阻力的抛射体运动模型，最后解得的是非线性、超越方程，因此如何解非线性方程是所需要知道的
逐步搜索法
二分法
简单迭代法(不动点迭代法)
设\(f(x)=0\)，\(f(x)\)在有根区间\([a,b]\)上是连续函数。
设计一个迭代公式，将\(f(x)=0\)写成等价形式:
在\([a,b]\)上任取\(x_0\)，带入上式右端，记所得的值为\(x_1\)
\[ x_1 = \varphi(x_0) \]同理得\(x_2 = \varphi(x_1), x_3 = \varphi(x_2)\)\
于是迭代公式为:
若
\[\lim_{n \to \infty} x_{n+1} = \lim_{n \to \infty} \varphi(x_n) = x^* = \varphi(x^*) \]称\(x_n\)为第\(n\)次近似
注：
- 当\(x_{n + 1}与x_n\)的差值的绝对值小于等于一个数时，认为是收敛了
- 迭代公式可能收敛,可能发散
- 迭代公式的设计方法可以有多种，但存在某种迭代公式是发散的，即无法求得解，因此迭代公式的设计方法要合适
收敛:
发散:
例：求\(e^{2x} + x - 4 = 0\)的根
将原方程改写为：
\[ e^{2x} = 4 - x \\ \Rightarrow 2x = \ln(4 - x) \\ \Rightarrow x = \cfrac{1}{2}\ln(4 - x) \]迭代公式即为：
\[ x_{n + 1} = \cfrac{1}{2}\ln(4 - x_n) \, , n = 0,1,2,3,... \]程序求解：
#include <iostream>
#include <cmath>
#include <vector>
#define EPS 1e-6
using namespace std;
int main(){
cout << "process:" << endl; //迭代过程：
vector<double> x;
x.push_back(0.0); // 等价于 x[0] = 0
cout << x[0] << endl;
int n;
for(n = 0; ; ++n){
x.push_back(0.5 * log(4 - x[n])); // 等价于 x[n + 1] = 0.5 * log(4 - x[n])
cout << x[n + 1] << endl;
if(abs(x[n + 1] - x[n]) <= EPS) // x[n + 1]与x[n]的差值小于等于EPS时，认为是收敛了
break;
}
cout << "times = " << n + 1 << endl;
cout << "result = " << x[n + 1] << endl; //输出结果
return 0;
}
结果：
process:
0
0.693147
0.597998
0.612182
0.610093
0.610401
0.610356
0.610362
0.610361
times = 8
result = 0.610361
另一种形式：
\[ x = 4 - e^{2x} \]迭代公式即为：
\[ x_{n + 1} = 4 - e^{2x_n} \, , n = 0,1,2,3,... \]程序求解(主要代码)：
将x.push_back(0.5 * log(4 - x[n]));改为x.push_back(4 - exp(2*x[n]));
结果：
process:
0
3
-399.429
4
-2976.96
4
-2976.96
4
-2976.96
4
-2976.96
4
-2976.96
4
-2976.96
4
-2976.96
4
...
程序陷入死循环，显然是发散，不收敛的，无法求得解
迭代函数收敛判定
若\(\varphi(x)\)在\(x^*\)的某个邻域内有一阶连续导数，且对该邻域内的\(x\)有\(|\varphi^\prime(x)|\leq q < 1\),则由微分中值定理得：
\[ |x_{n+1} - x^*| = |\varphi^\prime(\xi)(x_n - x^*)| \leq q|x_n - x^*| \]反复递推得：
\[ |x_n - x^*| \leq q|x_{n-1} - x^*| \leq q^2|x_{n-2} - x^*| \leq ... \leq q^n|x_0 - x^*| \]同时
\[ \begin{aligned} &\because q < 1 \\ &\therefore 当\, n \to \infty \, , q^n \to 0\\ &即 |x_n - x^*| \to 0 \, , x_n \to x^* \end{aligned} \]因此，递推公式若是收敛的，则需满足其在某一领域内有一阶连续导数，且导数的绝对值小于1
迭代函数收敛定理
综上，有以下定理：
设\(x^*\)是方程\(x=\varphi(x)\)的根，若\(\varphi(x)\)在\(x^*\)的某个邻域内有一阶导，且对该邻域内的一切\(x\)有
则迭代公式\(x_{n + 1} = \varphi(x_n)\)对该邻域内任一初值\(x_0\)均收敛
迭代常用条件 \(|x_n - x_{n+1}| < \varepsilon\)
在上面的例子中：
第一个迭代公式：
因此是收敛的
第二个迭代公式：
\[ x_{n + 1} = 4 - e^{2x_n} \\ |\varphi^\prime(x_n)| = 2e^{2x_n} \,不满足小于1 \,(x^*的领域内) \]因此是发散的
迭代法中有些迭代虽收敛，但速度慢!
加速算法 - 埃特金(Aitken)法
方程\(x = \varphi(x)\)如图所示：
根据初值\(x_0\)，得到\(A(x_0,\bar{x_1})，其中\bar{x_1} = \varphi(x_0)\)
再由\(\bar{x_1}\)，得到\(B(\bar{x_1},\bar{x_2})，其中\bar{x_2} = \varphi(\bar{x_1})\)
点\(A\)与\(B\)的连线与\(y = x\)得到一交点\(C(x_1,x_1)\)
根据\(x_1\)得到\((x_1,\bar{x_3})\)、\((\bar{x_3},\bar{x_4})\)、两者的连线与\(y=x\)的交点\((x_2,x_2)\)
根据\(x_2\) ... 得到\(x_3\)
根据\(x_3\) ... 得到\(x_4\)
...
即：
第一次，计算：
连接\(A(x_0,\bar{x_1}),B(\bar{x_1},\bar{x_2})\)，与\(y=x\)交于点\(C(x_1,x_1)\)，则有\(AC\)斜率 \(=\) \(CB\)斜率：
\[ \cfrac{x_1 - \bar{x_1}}{x_1 - x_0} = \cfrac{\bar{x_2} - \bar{x_1}}{\bar{x_1} - x_0} \\ \Rightarrow x_1 = \cfrac{x_0\bar{x_2} - \bar{x_1}^2}{x_0 - 2\bar{x_1} + \bar{x_2}} \]第二次，计算：
\[ \bar{x_1} = \varphi(x_1) \\ \bar{x_2} = \varphi(\bar{x_1}) \]同理可得：
\[ x_2 = \cfrac{x_1\bar{x_2} - \bar{x_1}^2}{x_1 - 2\bar{x_1} + \bar{x_2}} \]第三次：
\[ \bar{x_1} = \varphi(x_2) \\ \bar{x_2} = \varphi(\bar{x_1}) \\ x_3 = \cfrac{x_2\bar{x_2} - \bar{x_1}^2}{x_2 - 2\bar{x_1} + \bar{x_2}} \]...
综上，
第\(n + 1\)次，\(n = 0,1,2,...\)：
收敛定理
与函数迭代法一致
例：求\(e^{2x} + x - 4 = 0\)的根
将原方程改写为：
\[ x = \cfrac{1}{2}\ln(4 - x) \]程序求解(主要代码)：
vector<double> x;
x.push_back(0.0); // 等价于x[0] = 0
double bar_x1, bar_x2; // 对应\bar{x_1}和\bar{x_2}
cout << x[0] << endl;
int n;
for(n = 0; ; ++n){
bar_x1 = 0.5 * log(4 - x[n]); //计算\bar{x_1}
bar_x2 = 0.5 * log(4 - bar_x1); //计算\bar{x_2}
x.push_back((x[n] * bar_x2 - pow(bar_x1, 2))/(x[n] - 2 * bar_x1 + bar_x2)); //计算x[n + 1]
cout << x[n + 1] << endl;
if(abs(x[n + 1] - x[n]) <= EPS) // x[n + 1]与x[n]的差值小于等于EPS时，认为是收敛了
break;
}
结果：
process:
0
0.609483
0.610362
0.610362
times = 3
result = 0.610362
直观的可见，相比于普通函数迭代法(见上)，普通函数迭代次数为\(8\)，而经过加速算法，迭代次数缩减至\(3\)，收敛速度要快，精度也在\(1e-6\)
注：对某些发散过程，埃特金法也适用
牛顿迭代法(切线法)
设方程\(f(x)=0\)
首先在根\(x^*\)的附近取一点\(x_0\)作为初始近似值，过曲线\(y=f(x)\)上的点\((x_0,f(x_0))\)作切线
切线方程：
\[ y = f(x_0) + f^\prime(x_0)(x - x_0) \]第一次：
若\(f^\prime(x_0) \neq 0\)，该切线与\(x\)轴的交点为
\(x_1\)作为根\(x^*\)的第一次近似值
第二次：
以\(x_1\)过曲线\(y=f(x)\)上的点\((x_1,f(x_1))\)作切线，若\(f(x_1) \neq 0\)，该切线与\(x\)轴交点：
\(x_2\)作为根\(x^*\)的第二次近似值
...
综上，
第\(n + 1\)次，\(n = 0,1,2,...\)：
当近似值越趋近于方程的解\(x^*\)时，所得斜率会越来越大，即趋近方程的解\(x^*\)的迭代速度会越来越快
相比于埃特金(Aitken)法，埃特金(Aitken)法计算量大
收敛判断
牛顿法仍是迭代法，牛顿法迭代函数为：
\[ \varphi(x) = x - \cfrac{f(x)}{f^\prime(x)} \]若\(f^\prime(x)\)在根\(x\)的某个邻域内不为零，且\(f^{\prime\prime}(x)\)存在，
牛顿迭代法中，要求\(f^\prime(x_n) \neq 0\)，因此
$f^\prime(x)$在根$x$的某个邻域内不为零
判断迭代函数的收敛性需要对\varphi(x)求导，即求\(|\varphi^\prime(x)|\)，因此$f^{\prime\prime}(x)$存在
由迭代函数的收敛性判断定理可知，若：
\[ |\varphi^\prime(x)| = |\cfrac{f(x)f^{\prime\prime}(x)}{[f^{\prime}(x)]^2}|\leq q_1 < 1 \]\(q_1\)为常数，则牛顿迭代法收敛
保证牛顿迭代法收敛：
用逐步搜索法，找到方程的有根区间\([a,b]\)(尽量小)
使得在\([a,b]\)上，\(f^{\prime}(x),f^{\prime\prime}(x)\)都不变号，即\(f(x)\)在\([a,b]\)上保持严格单调(由\(f^{\prime}(x)\)确定)和凹向不变(由\(f^{\prime\prime}(x)\)确定)
有以下四种情况:
- 可以看出，用牛顿迭代法求得的序列\({x_n}\)均是单调地趋于\(x^*\)，故牛顿迭代法是收敛的。
- 凡满足关系式
\(x_0\)均可作为初始值
例如图(1)，(4)取\(x_0 = b\)，图(2)，(3)取\(x_0 = a\)
牛顿迭代法收敛定理
设函数\(f(x)\)在\([a,b]\)上存在二阶导数，且满足下列条件:
- \(f(a)f(b) < 0\)
- \(f^{\prime}(x_0)\)在\([a,b]\)上不为零
- \(f^{\prime\prime}(x_0)\)在\([a,b]\)上不变号
- 对\(x_0 \in [a,b]\)，有\(f(x_0)f^{\prime\prime}(x_0)>0\)
则牛顿法迭代序列\({x_n}\)收敛于方程\(f(x)=0\)在\((a,b)\)内的唯一根\(x^*\)
使用牛顿法，初始值\(x_0\)的选取很重要，若在\([a,b]\)上，\(f^{\prime}(x),f^{\prime\prime}(x)\)的符号不易判别，如何选取初值\(x_0\)？
如何选取初值\(x_0\)
由：
\[ x_1 = x_0 - \cfrac{f(x_0)}{f^\prime(x_0)} \]得：
\[ x_1 - x^* = (x_0 - x^*) - \cfrac{f(x_0)}{f^\prime(x_0)} \]设第\(n\)次的近似值与方程的根\(x^*\)的误差为\(e_n\)：
\[ e_n = x_n - x^* \]即\(e_0 = x_0 - x^*\),\(e_1 = x_1 - x^*\)
那么：
利用一阶及零阶泰勒公式泰勒公式：
\[ 0 = f(x^*) = f(x_0) + f^\prime(x_0)(x^* - x_0) + \cfrac{1}{2}f^{\prime\prime}(\xi)(x^* - x_0)^2 \\ 0 = f(x^*) = f(x_0) + f^\prime(\eta)(x^* - x_0) \]那么：
\[ \cfrac{e_1}{e_0} \approx \cfrac{-\cfrac{1}{2}f^{\prime\prime}(\xi)(x^* - x_0)^2}{f^\prime(x_0)(x^* - x_0)} = -\cfrac{1}{2}\cfrac{f^{\prime\prime}(\xi)(x^* - x_0)}{f^\prime(x_0)} = -\cfrac{1}{2}\cfrac{f^{\prime\prime}(\xi)-\cfrac{f(x_0)}{f^\prime(\eta)}}{f^\prime(x_0)} = \cfrac{f^{\prime\prime}(\xi)f(x_0)}{2f^\prime(x_0)f^\prime(\eta)} \]在\(x_0\)处可计算\(f(x),f^{\prime}(x_0),f^{\prime\prime}(x_0)\)的值，但无法计算\(f^{\prime\prime}(\xi)\)与\(f^\prime(\eta)\)
若\(f^{\prime}(x),f^{\prime\prime}(x)\)在\(x_0\)附近相对变化不大，即只要\(f^{\prime\prime}(x_0) \neq 0\)，则有近似公式:
要使牛顿迭代法收敛，误差必须减少，即\(|e_1|<|e_0|\)
\[ [f^\prime(x_0)]^2 > |\cfrac{f^{\prime\prime}(x_0)}{2}|\cdot|f(x_0)| \tag{3.2} \]只要\(x_0\)满足上式，且\(f^{\prime\prime}(x_0) \neq 0\)，就可将\(x\)作为牛顿迭代法初值
使用牛顿迭代法步骤
- 选初值\(x_0\)(为保证迭代收敛，可用收敛定理或式\((3.2)\)选定初值)，并计算\(f(x_0), f^{\prime}(x_0)\)
- 迭代，\(x_1 = x_0 - \cfrac{f(x_0)}{f^\prime(x_0)}\)，计算\(f(x_1), f^{\prime}(x_1)\)
- 若\(|x_1 - x_0| < \varepsilon\)(或\(|f(x_1)| < \varepsilon\))，迭代终止，\(x_1\)即为根，否则到
4. - 以\(x_1, f(x_1), f^{\prime}(x_1)\)替代\(x_0, f(x_0), f^{\prime}(x_0)\)回到
2.继续，直到\(|x_1 - x_0| < \varepsilon\)(或\(|f(x_1)| < \varepsilon\))为止
例：求方程\(f(x) = x^3 - 2x^2 + 4x + 1 = 0\)在\(x = 0\)附近的一个实根
对\(f(x)\)求导
\[ f^{\prime}(x) = 3x^2 - 4x + 4 \\ f^{\prime\prime}(x) = 6x - 4 \]程序求解(主要代码)：
vector<double> x;
x.push_back(1.0); // 初值x[0]
double x1, f_x, f_1_x; // 对应x_1, f(x_1), f^{\prime}(x_1)
cout << x[0] << endl;
int n;
for(n = 0; ; ++n){
f_x = pow(x[n], 3) - 2 * pow(x[n], 2) + 4 * x[n] + 1;
if(abs(f_x) <= EPS) //f(x)的绝对值小于等于EPS时，认为是收敛了
break;
f_1_x = 3 * pow(x[n], 2) - 4 * x[n] + 4;
x1 = x[n] - f_x/f_1_x;
x.push_back(x1);
cout << x[n + 1] << endl;
if(abs(x[n + 1] - x[n]) <= EPS) // x[n + 1]与x[n]的差值的绝对值小于等于EPS时，认为是收敛了
break;
}
结果：
process:
1
-0.333333
-0.228758
-0.222515
-0.222495
times = 5
result = -0.222495
初值取0.0结果：
process:
0
-0.25
-0.222892
-0.222495
times = 4
result = -0.222495
例：求\(e^{2x} + x - 4 = 0\)在区间\([0.5, 1]\)内的根的近似值
设：
\[ f(x) = e^{2x} + x - 4 \]则：
\[ f^{\prime}(x) = 2e^{2x} + 1 > 0 \\ f^{\prime\prime}(x) = 4e^{2x} > 0 \]且：
\[ f(1) = e^2 + 1 - 4 > 0 \]因此：取初值为\(x_0 = 1\)
程序求解(主要代码)：
修改f_x：f_x = exp(2 * x[n]) + x[n] - 4;
修改f_1_x：f_1_x = 2 * exp(2 * x[n]) + 1;
结果：
process:
1
0.721826
0.620693
0.610454
0.610362
times = 5
result = 0.610362
可见，使用牛顿迭代法也比一般的函数迭代法(\(8\)次)收敛速度要快
例：计算\(f(x) = x^{41} + x^3 + 1 = 0\)在\(x_0 = -1\)附近的实根
\[ f^{\prime}(x) = 41x^{40} + 3x^2 \\ \cfrac{1}{2}f^{\prime\prime}(x) = 820x^{39} + 3x \]可知：
\[ f(-1) = -1 \\ f^{\prime}(-1) = 44 \\ \frac{1}{2}f^{\prime\prime}(-1) = -823 \]\[ [f^{\prime}(-1)]^2 = 1936 > |\frac{f^{\prime\prime}(-1)}{2}|\cdot|f(-1)| = 823 \]因此，可取\(x_0 = -1\)为初值
...
例：计算平方根
\[ x = \sqrt{c} \]求\(f(x) = x^2 - c = 0\)的正根
\[ f^{\prime}(x) = 2x\\ f^{\prime\prime}(x) = 2 \\ x_{n + 1} = x_n - \cfrac{x_n^2 - c}{2x_n} = \cfrac{1}{2}(x_n + \cfrac{c}{x_n}) \]初值的选取：符合定理或条件
用于求复根
牛顿法可求方程的实根，也可以求方程的复根
令\(z_0 = x_0 + y_0i\)，迭代按：
例：计算\(f(z) = z^3 - z - 1 = 0\)的复根
取初值\(z_0 = -0.5 + 0.5i\)
...
弦截法(割线法)
牛顿法计算时要用到函数的导数，很多情况下难以使用，更何况选取初值还需使用到二阶导数
弦截法(割线法)是它的一种改进，在牛顿迭代公式中，用差商代替导数，即
\[ f^{\prime}(x_n) \approx \cfrac{f(x_n) - f(x_{n-1})}{x_n - x_{n-1}} \]代入牛顿迭代公式可得：
\[ x_{n+1} = x_n - \cfrac{f(x_n) }{f(x_n) - f(x_{n-1})}(x_n - x_{n - 1}) \, , n = 0,1,2,... \tag{4} \]迭代格式中没有用到函数的导数，计算方便
但收敛速度较牛顿迭代法要慢，开始时要用到两个不同的根的近似值作初值
例：求方程\(e^{2x} + x - 4 = 0\)在区间\([0.5,1]\)内根的近似值
设\(f(x) = e^{2x} + x - 4\)
两个初值取区间端点，\(x_0 = 0.5, x_1 = 1\)
使用弦截法步骤
- 选定初始值\(x_0, x_1\)，计算\(f(x_0),f(x_1)\)
- 按迭代公式 \(x_{n+1} = x_n - \cfrac{f(x_n) }{f(x_n) - f(x_{n-1})}(x_n - x_{n - 1}) \, , n = 0,1,2,...\)
计算\(x_2\)，再计算\(f(x_2)\) - 判定若\(|f(x_2)| < \varepsilon\)，则迭代终止
否则，用\((x_2, f(x_2))和(x_1, f(x_1)),\)分别代替\((x_1, f(x_1)),(x_0, f(x_0))\)
重复步骤2，3，直至\(|f(x_2)| < \varepsilon\)