【数值计算方法】《工程数值计算Python教程》笔记
第一章:绪论
1.1 1.1 1.1|数值计算在工程科学中的重要性
1.2 1.2 1.2|数值计算方法
1.3 1.3 1.3|程序设计
盒图
- 盒图也称为 N ? S N-S N?S结构化流程图
- 如果一个算法可以用 N ? S N-S N?S流程图描述,则说明它是结构化的
- 盒图示例
计算方法的选取
减少运算次数
- 计算多项式
P
n
(
x
)
=
a
0
+
a
1
x
+
a
2
x
2
+
?
+
a
n
x
n
P_{n} (x) = a_{0} + a_{1} x + a_{2} x^{2} + \cdots + a_{n} x^{n}
Pn?(x)=a0?+a1?x+a2?x2+?+an?xn的值
- 如果逐项计算,所需乘法次数为 1 + 2 + ? + n = 1 2 n ( n + 1 ) 1 + 2 + \cdots + n = \cfrac{1}{2} n (n + 1) 1+2+?+n=21?n(n+1),所需加法次数为 n n n
S 0 = a 0 S k = a k x k ( k = 1 , 2 , ? ? , n ) P n ( x ) = ∑ i = 0 n S i S_{0} = a_{0} \\ S_{k} = a_{k} x^{k} (k = 1 , 2 , \cdots , n) \\ P_{n} (x) = \displaystyle\sum\limits_{i = 0}^{n}{S_{i}} S0?=a0?Sk?=ak?xk(k=1,2,?,n)Pn?(x)=i=0∑n?Si?
-
- 如果采用递推法(秦九韶算法),则只需 n n n次乘法和 n n n次加法
S n = a n S k = x S k ? 1 + a k ( k = n ? 1 , n ? 2 , ? ? , 0 ) P n ( x ) = S 0 S_{n} = a_{n} \\ S_{k} = x S_{k - 1} + a_{k} (k = n - 1 , n - 2 , \cdots , 0) \\ P_{n} (x) = S_{0} Sn?=an?Sk?=xSk?1?+ak?(k=n?1,n?2,?,0)Pn?(x)=S0?
避免相近的数相减
- 在计算过程中,应避免相近的数相减,不然会使有效数字的位数大大减少
1.4 1.4 1.4|误差的来源、表示及传递
误差的来源和分类
模型误差
- 模型与实际问题的误差
观测误差
- 通常数学模型中都包含一些需要实验测定的参数,这些参数存在测定值与真实值之间的误差
截断误差
- 数值计算方法中,常用收敛无穷级数的前几项代替无穷级数
- 截断误差与算法有关,常用截断误差限或截断误差的阶来判断某种算法的优劣
舍入误差
- 由于计算机位数有限引起的误差称为舍入误差
- 舍入误差与计算机有关,也与数学表达式有关
误差的表示
绝对误差
- 用 x x x表示准确值 x ? x^{*} x?的近似值,绝对误差为 e = x ? x ? e = x - x^{*} e=x?x?
- 通常 x ? x^{*} x?是未知的,故 e e e的真实值也是不知道的,常根据实际情况估算它的上限, ∣ e ∣ = ∣ x ? x ? ∣ ≤ ε | e | = | x - x^{*} | \leq \varepsilon ∣e∣=∣x?x?∣≤ε, x ? ε ≤ x ? ≤ x + ε x - \varepsilon \leq x^{*} \leq x + \varepsilon x?ε≤x?≤x+ε, ε \varepsilon ε称为 x x x的绝对误差限
相对误差
- e r = e x = x ? x ? x e_{r} = \cfrac{e}{x} = \cfrac{x - x^{*}}{x} er?=xe?=xx?x??
- ∣ e r ∣ = ∣ x ? x ? x ∣ ≤ ε r | e_{r} | = | \cfrac{x - x^{*}}{x} | \leq \varepsilon_{r} ∣er?∣=∣xx?x??∣≤εr?
平均误差
- 算数平均值 x ˉ = 1 n ∑ i = 1 n x i \bar{x} = \cfrac{1}{n} \displaystyle\sum\limits_{i = 1}^{n}{x_{i}} xˉ=n1?i=1∑n?xi?
- 平均误差为 δ = 1 n ∑ i = 1 n d i = 1 n ∑ i = 1 n ∣ x ˉ ? x i ∣ \delta = \cfrac{1}{n} \displaystyle\sum\limits_{i = 1}^{n}{d_{i}} = \cfrac{1}{n} \displaystyle\sum\limits_{i = 1}^{n}{| \bar{x} - x_{i} |} δ=n1?i=1∑n?di?=n1?i=1∑n?∣xˉ?xi?∣
标准误差
- σ = 1 n ? 1 ∑ i = 1 n d i 2 \sigma = \sqrt{\cfrac{1}{n - 1} \displaystyle\sum\limits_{i = 1}^{n}{d_{i}^{2}}} σ=n?11?i=1∑n?di2??
- 标准误差也称为均方根误差
误差的传递
误差在和、差计算中的传递
- 设 x ? x^{*} x?、 y ? y^{*} y?的近似值为 x x x、 y y y, z ? = x ? + y ? z^{*} = x^{*} + y^{*} z?=x?+y?
绝对误差
- e z = e x + e y e_{z} = e_{x} + e_{y} ez?=ex?+ey?
- ∣ e z ∣ ≤ ∣ e x ∣ + ∣ e y ∣ | e_{z} | \leq | e_{x} | + | e_{y} | ∣ez?∣≤∣ex?∣+∣ey?∣
相对误差
-
e
r
z
=
(
x
+
y
)
?
(
x
?
+
y
?
)
x
+
y
=
x
?
x
?
x
×
x
x
+
y
+
y
?
y
?
y
×
y
x
+
y
e_{rz} = \cfrac{(x + y) - (x^{*} + y^{*})}{x + y} = \cfrac{x - x^{*}}{x} \times \cfrac{x}{x + y} + \cfrac{y - y^{*}}{y} \times \cfrac{y}{x + y}
erz?=x+y(x+y)?(x?+y?)?=xx?x??×x+yx?+yy?y??×x+yy?
- 当 x x x与 y y y同号时, ∣ x x + y ∣ ≤ 1 | \cfrac{x}{x + y} | \leq 1 ∣x+yx?∣≤1, ∣ y x + y ∣ ≤ 1 | \cfrac{y}{x + y} | \leq 1 ∣x+yy?∣≤1, e r z ≤ ∣ e r x ∣ + ∣ e r y ∣ e_{rz} \leq | e_{rx} | + | e_{ry} | erz?≤∣erx?∣+∣ery?∣
- 当 x + y ≈ 0 x + y \approx 0 x+y≈0时,可能有 e r z ? ∣ e r x ∣ + ∣ e r y ∣ e_{rz} \gg | e_{rx} | + | e_{ry} | erz??∣erx?∣+∣ery?∣
误差在积、商计算中的传递
- 设 x ? x^{*} x?, y ? y^{*} y?均为正数,近似值分别为 x x x、 y y y,绝对误差为 d x = x ? x ? dx = x - x^{*} dx=x?x?, d y = y ? y ? dy = y - y^{*} dy=y?y?,相对误差为 x ? x ? x = d x x = d ln ? x \cfrac{x - x^{*}}{x} = \cfrac{dx}{x} = d \ln{x} xx?x??=xdx?=dlnx, y ? y ? y = d y y = d ln ? y \cfrac{y - y^{*}}{y} = \cfrac{dy}{y} = d \ln{y} yy?y??=ydy?=dlny
乘积的绝对误差
- d ( x y ) = x d y + y d x d(xy) = x dy + y dx d(xy)=xdy+ydx
乘积的相对误差
- d ln ? x y = d ( ln ? x + ln ? y ) = d x x + d y y d \ln{xy} = d(\ln{x} + \ln{y}) = \cfrac{dx}{x} + \cfrac{dy}{y} dlnxy=d(lnx+lny)=xdx?+ydy?
商的绝对误差
- d ( x / y ) = y d x ? x d y y 2 d(x / y) = \cfrac{y dx - x dy}{y^{2}} d(x/y)=y2ydx?xdy?
商的相对误差
- d ln ? x / y = d ( ln ? x ? ln ? y ) = d x x ? d y y d \ln{x / y} = d(\ln{x} - \ln{y}) = \cfrac{dx}{x} - \cfrac{dy}{y} dlnx/y=d(lnx?lny)=xdx??ydy?
第二章:Python
基础
2.1 2.1 2.1|概述
第三章:方程(组)的求解
3.1 3.1 3.1|非线性代数方程的求根
二分法
原理
- 设 f f f是一个连续函数,假设存在一个区间 [ a , b ] [a,b] [a,b], f f f在区间的两端点处的值符号相反,即 f ( a ) f ( b ) < 0 f(a) f(b) < 0 f(a)f(b)<0,于是, f f f在区间 [ a , b ] [a,b] [a,b]有一个根,这可由函数的中值定理得到
步骤
-
有一个区间 [ a , b ] [a,b] [a,b]和值 u = f ( a ) u = f(a) u=f(a)、 v = f ( b ) v = f(b) v=f(b),满足 u v < 0 uv < 0 uv<0
-
构造区间的中点 c = ( a + b ) / 2 c = (a + b) / 2 c=(a+b)/2并计算 w = f ( c ) w = f(c) w=f(c),可能碰巧 w = 0 w = 0 w=0,则算法结束
-
通常情况下 w ≠ 0 w \neq 0 w=0,如果 w u < 0 wu < 0 wu<0,则可断定在区间 [ a , c ] [a,c] [a,c]中存在 f f f的一个根,因而将 c c c的值存储到 b b b中,将 w w w的值存储到 v v v中,如果 w v < 0 wv < 0 wv<0,则 f f f在区间 [ c , b ] [c,b] [c,b]中有一个根,因而将 c c c的值存储到 a a a中,将 w w w的值存储到 u u u中
-
重复进行,每次区间长度减少一半,直到区间达到需要的精度
-
结束时,根的最佳估计是 ( a + b ) / 2 (a + b) / 2 (a+b)/2
-
在程序运行过程中,中点处函数值 w = 0 w = 0 w=0的概率是极低的,为此增加一个判断语句是不值得的,对任意区间 [ a , b ] [a,b] [a,b],只要 u v ≤ 0 uv \leq 0 uv≤0即可保证区间中包含零点
Python
实现
- 求解 f ( x ) = x 3 ? 3 x + 1 f(x) = x^{3} - 3x + 1 f(x)=x3?3x+1在区间 [ 0 , 1 ] [0,1] [0,1]上的零点
def bisect(f, a, b, eps=1e-9, args=()):
u, v = f(a, *args), f(b, *args)
if u * v <= 0:
while abs(a - b) > eps:
c = (a + b) / 2
w = f(c, *args)
if w * u <= 0:
b, v = c, w
else:
a, u = c, w
return (a + b) / 2
else:
print(f'The function values at {a} and {b} have the same sign')
def f(x): return x * (x * x - 3) + 1
if __name__ == '__main__':
a, b = 0, 1
eps = 1e-9
res = bisect(f, a, b, eps)
print('The result is:', res)
- 求解
f
(
x
)
=
cosh
?
(
x
2
+
1
?
e
x
)
+
log
?
∣
sin
?
x
∣
f(x) = \cosh{(\sqrt{x^{2} + 1} - e^{x})} + \log{| \sin{x} |}
f(x)=cosh(x2+1??ex)+log∣sinx∣在区间
[
0
,
1
]
[0,1]
[0,1]上的零点
- 该式在 0 0 0点没有定义,计算会发生溢出,此时可以选择一个很小的值作为左边界,如 1 0 ? 10 10^{-10} 10?10
import numpy as np
from bisection import bisect
def f(x): return np.cosh(np.sqrt(x * x + 1) - np.exp(x)) + np.log10(np.abs(np.sin(x)))
if __name__ == '__main__':
a, b = 1e-10, 1
res = bisect(f, a, b)
print(f'The result is: {res}')
本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。 如若内容造成侵权/违法违规/事实不符,请联系我的编程经验分享网邮箱:veading@qq.com进行投诉反馈,一经查实,立即删除!