【数值计算方法】一、非线性方程的数值解法（二分法+迭代法+Python实现）

大家好，欢迎来到IT知识分享网。

立个 flag 在这里，暑假学完数值分析方法（本科版），等我学完后再回这 flag 来看看。

[参考教材]：《数值计算方法》马东升董宁

[参考视频]：数值计算方法数值分析计算方法_哔哩哔哩_bilibili

一、非线性方程组定义

方程 f(x)=0 当 f(x)是一次多项式时，称f(x)为线性方程。

代数多项式方程

$f(x)=a_{n}x^{n}+a_{n-1}x^{n-1}+\cdot\cdot\cdot+a_{1}x+a_{0}=0$

超越方程

$e^{-x}-sin(\frac{\pi x}{2})=0$

重根

若 $f(x)=(x-x^{*})^{m}g(x)=0,g(x)\neq 0$ ，m为正整数，则称 $x^{*}$ 为 f(x)=0 的 m 重根。若函数 f(x) 有 m 阶连续导数，方程 f(x)=0 的充要条件是：

$f(x^{*})=f^{'}(x^*)=\cdot\cdot\cdot=f^{m-1}(x^*)=0,f^m\neq 0$

二、方程求根特点

1、根的存在性

n次代数方程有n个根

2、根的分布

在 [a,b] 内连续，且 $f(a)\times f(b)<0$ （即f(a）和f(b)异号)，则在[a, b]内 f(x)=0 至少有一个根。（有根区间）
在 [a,b] 内连续，严格单调，且 $f(a)\times f(b)<0$ （即f(a）和f(b)异号)，则在[a, b]内 f(x)=0 有且仅有一个根。（隔根区间）

确定隔根区间的方法：

作 y=f(x) 的草图，看 f(x) 与 x 轴交点定 [a, b]
逐步搜索，在连续区间 [a, b] 内，选取适当的 $x_1,x_2\subseteq (a,b)$ ，若 $f(a)\times f(b)<0$ ，则内有根。
根的精确化，将区间逐步缩小，直到找到更加精确的根。

三、二分法

1、区间对分

将区间对分，判别 f(x) 的符号，逐步缩小有根区间。

2、方法

流程图如下，其中 $\varepsilon$ 为预先给定的精度

重复上述流程，直到找到方程的根。

3、收敛性

由二分法产生一个有根区间:

$[a,b]\supset [a_1,a_2]\supset[a_2,b_2]\supset\cdot\cdot\cdot\supset[a_k,b_k]$

关于的区间长度：

$b_k-a_k=\frac{1}{2}(b_{k-1}-a_{k-1})=\cdot\cdot\cdot=\frac{1}{2^k}(b-a)$

当 k 足够大时，取近似值 $x_k=\frac{a_k+b_k}{2}$ 为方程的根

误差： $\left |x^*- x_k \right |\leqslant \frac{b_k-a_k}{2}=\frac{b-a}{2^{k+1}}<\varepsilon$

先验估计，设 0″> 为给定精度要求，可确定分半次数k使得 $\left |x^*- x_k \right |\leqslant \frac{b-a}{2^{k+1}}<\varepsilon$

两边取对数得到 \frac{ln(b-a)-ln(\varepsilon ))}{ln2}-1″ />

4、例题

用二分法求在 [1, 2] 内的根，要求绝对误差不超过 $\frac{1}{2}\times10^{-2}$

求出，f(1)=-5, f(2)=14

k	有根区间	中点	f(x)
0	[1, 2]	1.5	f(1.5)=2.375>0
1	[1, 1.5]	1.25	f(1.25)=-1.<0
2	[1.25, 1.5]	1.375	f(1.375)=0.>0
3	[1.25, 1.375]	1.313	f(1.313)=-0.<0
4	[1.313, 1.375]	1.344	f(1.344)=-0.<0
5	[1.344, 1.375]	1.360	f(1.360)=-0.086144<0
6	[1.360, 1.375]	1.368	f(1.368)=0.045804>0
7	[1.360, 1.368]	1.364	f(1.364)=-0.020299<0

上述表格中： $\frac{b_k-a_k}{2}<\varepsilon =\frac{1}{2}\times10^{-2}=0.005$

$\frac{1.368-1.360}{2}=0.004<0.005$

ps:要求精确到小数后第三位，即使要求 $\left | x-x_k \right |<\frac{1}{2}\times10^{-3}$

5、二分法评价

优点：计算简单，方法可靠，只要求f(x)连续
缺点：不能求重根，不能求复根，收敛速度不算太快

在求方程近似根时，一般不单独使用，常用来为其他方法提供初值

6、代码实现

挺简单的，改一下fx就好

""" @ S_Iris method of bisection """ def fx(x): ans = x3+4*(x2)-10 return ans def bisection(fx, a, b, epsilon): if a < b: ak = a bk = b elif a == b: return a else: ak = b bk = a k = 0 xmid = ak while (bk - ak)/2 >= epsilon: xmid = ak + (bk - ak)/2 ans = fx(xmid) if ans < 0: ak = xmid elif ans == 0: return xmid else: bk = xmid xmid = ak + (bk - ak) / 2 return xmid def main(): a = 1 b = 2 epsilon = 0.005 x = bisection(fx, a, b, epsilon) print("方程的近似解为{}".format(x)) if __name__ == "__main__": main()

四、迭代法

1、迭代法的一般过程

改写方程 f(x)=0 为 x=g(x)形式，要求g(x)连续且收敛
建立迭代格式： $x_{n+1}=g(x_n)$ ，得到序列 $\left \{ x_n \right \}$
若 $\left \{ x_n \right \}$ 收敛，必收敛到 f(x)=0 的根：
$\lim_{n \to \infty }x_{n+1}=\lim_{n \to \infty }g(x_{n})=g(\lim_{n \to \infty}x_n)$
若 $\left \{ x_n \right \}$ 收敛，即 $\lim_{n\to\infty}x_n=x^*$ ，则：
$x^*=g(x^*)\Rightarrow f(x^*)=0$

2、几何表示

$x=g(x)\Rightarrow \left\{\begin{matrix}y=g(x) \\x=y \end{matrix}\right.$ 交集即为真根

3、例题

【例】用迭代法求方程的根

解法一：

方程改写为 $x=\sqrt{2x+3}$
建立迭代公式 $x_{k+1}=\sqrt{2x_k+3},(k=0,1,2,3,\cdot\cdot\cdot)$
取（设初值），则根据迭代公式有
当 $k \to \infty$ 时， $x_k\to3$

解法二（反例）：

方程改写成 $x=\frac{1}{2}(x^2-3)$
建立迭代公式 $x_{k+1}=\frac{1}{2}(x_k^2-3),(k=0,1,2,3,\cdot\cdot\cdot)$
取（设初值），则根据迭代公式有，当 $k\to\infty$ ， $x_k\to\infty$ ，发散。

从以上两个例子可以看出，迭代法必须要求迭代函数的导数满足条件： $\left | g'(x) \right |<1$

4、整体收敛性

考虑方程 x=g(x)，若：

当 $x\in [a, b]$ 时， $g(x)\in[a, b]$ ；
${\exists}L(0<L<1)$ 使得 $\left | g'(x) \right |\leqslant L<1$ 对 ${\forall}x\in[a,b]$ 成立。则任取 $x_0\in[a,b]$ ，由 $x_{k+1}=g(x_k)$ 得到的序列 $\left \{ x_k \right \}_{k=0}^{\infty}$ 收敛于 g(x) 在 [a, b]上的唯一不动点，并且有误差估计式：
（1） $\left | x^*-x_k \right |\leqslant \frac{1}{1-L}\left | x_{k+1}-x_k \right |$
（2） $\left | x^*-x_k \right |\leqslant \frac{L^*}{1-L}\left | x_{1}-x_0 \right |$
且存在极限 $\lim _{k\to\infty}\frac{x^*-x_{k+1}}{x^*-x_k}=g'(x^*)$

5、局部收敛性

如果函数 g(x) 在 x^* 的一个领域 $[x^*-\delta ^*, x^*+\delta ^*]$ 内连续可微， x^* 为方程 x=g(x) 的根，且 $\left | g'(x) \right |<1$ ，则存在正整数 $\delta,\delta \leqslant \delta ^*$ ，使得对任意 $x_0\in[x^*-\delta ,x^*+\delta ]$ ，迭代序列 $x_{k+1}=g(x_k),(k=0,1,2,3,...)$ 收敛于 x^*

6、代码实现

使用了 SymPy 库来进行收敛性判断，使用时需要注意数据格式

""" @ S_Iris iteration_method """ import sympy def iteration(expr, x_, x0, precision_): """ :param expr: gx表达式 :param x_: 变量符号 :param x0: 初值 :param precision_:精度要求 :return: 是否收敛（bool），根x* """ a = str(x_) precision = precision_ x = sympy.symbols(a) gx = expr print("gx={}".format(gx)) d_gx = sympy.diff(gx, x) x_k = x0 while abs((gx.subs(x, x_k)).evalf()-x_k) > precision: # 判断收敛性 if d_gx.subs(x, x_k).evalf() > 1: print("gx不收敛！！！") return False, 0 x_k = gx.subs(x, x_k) x_k = x_k.evalf() print("x*:{}".format(x_k)) return True, x_k def main(): fx = "x2-2*x-3" fx = sympy.sympify(fx) print("fx={}".format(fx)) x = sympy.symbols('x') expr1 = sympy.sqrt(2*x+3) expr2 = 0.5*((x2)-3) precision = 10(-9) print("gx表达式1") flag1, ans1 = iteration(expr1, x, 4.0, precision) print("gx表达式2") flag2, ans2 = iteration(expr2, x, 4.0, precision) if __name__ == "__main__": main()

运行结果：

可以看到与理论计算结果一致。

五、牛顿迭代法

1、原理

在迭代法中构造 g(x) 的一条重要途径：用近似方程代替原方程求根
牛顿法的思想是将非线性方程线性化。

2、方法(Taylor展开)

设 x_k 是 f(x)=0 的一个近似根，将 f(x) 在 x_k 出做一阶泰勒展开得到

$f(x)=f(x_k)+f'(x_k)(x-x_k)+\frac{f''(x_k)}{2!}(x-x_k)^2+...+\frac{f^{(n)}(x_k)}{n!}(x-x_k)^n$

舍掉余项得： $f(x)\approx f(x_k)+f'(x_k)(x-x_k)$

则有 f(x)=0 近似为

解出 $x = x_k-\frac{f(x_k)}{f'(x_k)}$

将写成 $x_{k+1}$ ，得到迭代公式

$x_{k+1} = x_k-\frac{f(x_k)}{f'(x_k)}$

3、几何意义

4、例题

【例】用牛顿法求 $f(x)=e^{\frac{-x}{4}}(2-x)-1=0$ 的根

解：

显然， $f(0)\times f(2)<0$ ，方程于 [0, 2] 内有一根。
求导 $f'(x)=\frac{e^{\frac{-x}{4}}(x-6)}{4}$
迭代公式 $x_{k+1}=x_k-\frac{e^{\frac{-x_k}{4}}(2-x_k)-1}{e^{\frac{-x_k}{4}}(x_k-6)/4},(k=0,1,2,...)$ ()

（1）取 x_0=1.0

k
0	1.0
1	-1.15599
2	0.
3	0.
4	0.
5	0.
6	0.

此时迭代公式收敛

（2）取 x_0=8.0

k
0	8.0
1	34.
2	869.1519

此时迭代公式发散

当初值选取靠近根的时候，牛顿法收敛快，当初值不是选取接近方程根时，牛顿法可能会给出发散的结果。

5、收敛定理（牛顿法）

收敛的充分条件，设 $f\in C^2[a, b]$ 若

; （有根）
在整个 [a, b] 上不变号且 $f'(x)\neq 0$ ;（根唯一）
选取 $x_0\in [a,b]$ 使得 0″> ;（保证收敛）

则用牛顿法产生的序列 $\left \{ x_k \right \}$ 收敛到 f(x) 在 [a, b] 的唯一根。

局部收敛：

$f\in C^2[a, b]$ ，若 x^* 为 f(x) 在 [a, b] 上的根，且 $f'(x^*)\neq 0$ ，则在 x^* 的邻域 $S\delta (x^*)$ 使得任取初值 $x_0\in S\delta (x^*)$ ，牛顿法产生的序列 $\left \{ x_k \right \}$ 收敛到 x^* ，且满足：

$\lim_{k\to \infty}\frac{x^*-x_{k+1}}{(x^*-x_k)^2}=-\frac{f''(x^*)}{2f'(x^*)}$

6、代码实现

由于牛顿法的收敛是充分条件，不是必要条件，不方便直接判断收敛性，所以我们设置一个最大迭代次数，如果超过最大迭代次数还未得到结果，则认为该初值条件发散。

""" @ S_Iris newton_iteration_method """ import sympy def newton_teration(expr, x_, x0, precision_, m): """ :param expr: gx表达式 :param x_: 变量符号 :param x0: 初值 :param precision_:精度要求 :param m: 最高迭代次数 :return: 是否收敛（bool），根x* """ a = str(x_) precision = precision_ x = sympy.symbols(a) fx = expr d_fx = sympy.diff(fx, x) d_2_fx = sympy.diff(d_fx, x) x_k = x0 x_k_add_1 = fx/d_fx n = 0 while abs((fx.subs(x, x_k)).evalf()) > precision and n < m: n += 1 print("\r迭代{}轮：".format(n), end="") x_k = x_k - x_k_add_1.subs(x, x_k).evalf() # print("第{}轮x_k值为：{}".format(n, x_k)) if n == m: print("未找到解，初值条件可能不收敛") return False, 0 print("找到方程的根") print("x*:{}".format(x_k)) return True, x_k def main(): x = sympy.symbols('x') fx = "exp(-x/4)*(2-x)-1" fx = sympy.sympify(fx) print("fx={}".format(fx)) precision = 10(-9) print("初值为0.0") flag1, ans1 = newton_teration(fx, x, 0.0, precision, 100) print("初值为1.0") flag1, ans1 = newton_teration(fx, x, 1.0, precision, 100) print("初值为2.0") flag1, ans1 = newton_teration(fx, x, 2.0, precision, 100) print("初值为3.0") flag1, ans1 = newton_teration(fx, x, 3.0, precision, 100) print("初值为5.9") flag1, ans1 = newton_teration(fx, x, 5.9, precision, 100) print("初值为5.95") flag1, ans1 = newton_teration(fx, x, 5.95, precision, 100) print("初值为5.") flag1, ans1 = newton_teration(fx, x, 5., precision, 100) # print("初值为6.0") # flag1, ans1 = newton_teration(fx, x, 6.0, precision, 100) # print("初值为6.0000001") # flag1, ans1 = newton_teration(fx, x, 6.0000001, precision, 100) # print("初值为6.001") # flag1, ans1 = newton_teration(fx, x, 6.001, precision, 100) # print("初值为8.0") # flag1, ans1 = newton_teration(fx, x, 8.0, precision, 100) if __name__ == "__main__": main()

运行结果，这个运行结果十分有趣，我尝试了几个初值条件，最后发现这个初值条件十分特殊。

直接看运行结果

我们可以看到在初值条件5.9及一下时，都能在100轮一下找到方程的解，而5.95则没有找到方程的解，但是当我将5.95的结果值打印出来后，初值5.95是收敛的，直接看图。

前20轮，中间省略了，太长了，直接展示后20轮

也就是说，当我的m设置得再大一点，初值5.95用牛顿迭代法就可能找到根

事实证明，在第177轮初值5.95找到了方程的根。

但遗憾的是，初值5.，在经历轮后依旧没有找到方程的根

而在初值条件为6.0时，程序报错

也就是程序在运算的时候出现了 NAN，也就是空，原因在于在牛顿法进行迭代的时候，我们需要一阶导数作为迭代公式分母，但是在 x_0=6 正好是一阶导数零点，所以不能作为分母进行运算。

在初值条件为6.0000001时，程序会卡死。

在初值条件为 6.1 的时候，程序报错

溢出错误，代表该数过大，无法继续参与云算，此时函数才算是正真意义上的发散。

在初值条件大于 6 时，其他数同理，可以试试。

免责声明：本站所有文章内容,图片，视频等均是来源于用户投稿和互联网及文摘转载整编而成，不代表本站观点，不承担相关法律责任。其著作权各归其原作者或其出版社所有。如发现本站有涉嫌抄袭侵权/违法违规的内容,侵犯到您的权益，请在线联系站长,一经查实,本站将立刻删除。本文来自网络,若有侵权，请联系删除，如若转载，请注明出处：https://haidsoft.com/126630.html

【数值计算方法】一、非线性方程的数值解法（二分法+迭代法+Python实现）

一、非线性方程组定义

代数多项式方程

超越方程

重根

二、方程求根特点

1、根的存在性

2、根的分布

三、二分法

1、区间对分

2、方法

3、收敛性

4、例题

5、二分法评价

6、代码实现

四、迭代法

1、迭代法的一般过程

2、几何表示

3、例题

4、整体收敛性

5、局部收敛性

6、代码实现

五、牛顿迭代法

1、原理

2、方法(Taylor展开)

3、几何意义

4、例题

5、收敛定理（牛顿法）

6、代码实现

相关推荐

发表回复