非线性方程的解
非线性方程的解
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_2 = \varphi(x_1), x_3 = \varphi(x_2)\)\
于是迭代公式为:
若
称\(x_n\)为第\(n\)次近似
注:
- 当\(x_{n + 1}与x_n\)的差值的绝对值小于等于一个数时,认为是收敛了
- 迭代公式可能收敛,可能发散
- 迭代公式的设计方法可以有多种,但存在某种迭代公式是发散的,即无法求得解,因此迭代公式的设计方法要合适
收敛:
发散:
例:求\(e^{2x} + x - 4 = 0\)的根
将原方程改写为:
迭代公式即为:
程序求解:
#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.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\),则由微分中值定理得:
反复递推得:
同时
因此,递推公式若是收敛的,则需满足其在某一领域内有一阶连续导数,且导数的绝对值小于1
迭代函数收敛定理
综上,有以下定理:
设\(x^*\)是方程\(x=\varphi(x)\)的根,若\(\varphi(x)\)在\(x^*\)的某个邻域内有一阶导,且对该邻域内的一切\(x\)有
则迭代公式\(x_{n + 1} = \varphi(x_n)\)对该邻域内任一初值\(x_0\)均收敛
迭代常用条件 \(|x_n - x_{n+1}| < \varepsilon\)
在上面的例子中:
第一个迭代公式:
因此是收敛的
第二个迭代公式:
因此是发散的
迭代法中有些迭代虽收敛,但速度慢!
加速算法 - 埃特金(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\)斜率:
第二次,计算:
同理可得:
第三次:
...
综上,
第\(n + 1\)次,\(n = 0,1,2,...\):
收敛定理
与函数迭代法一致
例:求\(e^{2x} + x - 4 = 0\)的根
将原方程改写为:
程序求解(主要代码):
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))\)作切线
切线方程:
第一次:
若\(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)法计算量大
收敛判断
牛顿法仍是迭代法,牛顿法迭代函数为:
若\(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)$存在
由迭代函数的收敛性判断定理可知,若:
\(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\)
由:
得:
设第\(n\)次的近似值与方程的根\(x^*\)的误差为\(e_n\):
即\(e_0 = x_0 - x^*\),\(e_1 = x_1 - x^*\)
那么:
利用一阶及零阶泰勒公式泰勒公式:
那么:
在\(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|\)
只要\(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)\)求导
程序求解(主要代码):
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]\)内的根的近似值
设:
则:
且:
因此:取初值为\(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\)附近的实根
可知:
因此,可取\(x_0 = -1\)为初值
...
例:计算平方根
求\(f(x) = x^2 - c = 0\)的正根
初值的选取:符合定理或条件
用于求复根
牛顿法可求方程的实根,也可以求方程的复根
令\(z_0 = x_0 + y_0i\),迭代按:
例:计算\(f(z) = z^3 - z - 1 = 0\)的复根
取初值\(z_0 = -0.5 + 0.5i\)
...
弦截法(割线法)
牛顿法计算时要用到函数的导数,很多情况下难以使用,更何况选取初值还需使用到二阶导数
弦截法(割线法)是它的一种改进,在牛顿迭代公式中,用差商代替导数,即
代入牛顿迭代公式可得:
迭代格式中没有用到函数的导数,计算方便
但收敛速度较牛顿迭代法要慢,开始时要用到两个不同的根的近似值作初值
例:求方程\(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\)