数值分析 笔记（一）数值计算中的误差

数值分析

专栏目录

---

0. 前言

数值分析 也称 计算数学， 为数学科学的一个分支。它主要对计算机求解各种数学问题的数值计算方法与实现。在实际应用中，根据问题建立的 数学模型 其 完备形式 一般不能方便得求出 精确解，对复杂的 非线性模型 忽略一些因素 简化模型 往往不能满足 近似程度的要求，因此使用 数值方法 直接求解 做 较少的简化，可以得到满足一定近似程度要求的结果。

import numpy as np
import re

1. 数值计算 的 误差

由于我们 对实际问题 进行抽象描述 简化 后建立数学模型，其与实际问题之间 是近似的，我们对他们之间出现的这种误差 称为 模型误差。由于 这种误差 难于用数量表示，假设 数学模型合理，我们在此不讨论此误差。在实际应用中 还存在 观测误差，其不为计算过程中生成的误差，同样在此不做讨论。

在数学模型不能获得精确解时，通过数值方法求得的近似解与精确解之间所导致的误差 为在此讨论的误差。由于计算机上的字长有限，原始数据在计算机上表示时就有误差，计算过程中又产生新的误差，这种误差称为 舍入误差。例：用 $3.14159$ 来近似代替 $\pi$, 产生的误差为

$$R=\pi -3.14159=0.0000026\cdots$$

1.1. 误差

• 绝对误差：$e^{\ast }=x^{\ast }-x$
  • $x$ 为 准确值
  • $x^{\ast }$ 为 近似值
• 误差限: 误差绝对值的一个上限，该值总为正数
  • ${\epsilon }^{\ast }$ 为 近似值的误差限 $|x^{\ast }-x|\le {\epsilon }^{\ast }$
  • 即 $x^{\ast }-{\epsilon }^{\ast }\le x\le x^{\ast }+{\epsilon }^{\ast }$
  • 有时 使用 $x=x^{\ast }\pm {\epsilon }^{\ast }$ 表示

def 误差(x_star: float, x_true:float):
    return x_star - x_true

• 相对误差:
  • 考虑误差大小时，其值对 $x$ 本身大小相关，因此应考虑它们的比值。
  • 近似值 $x^{\ast }$ 的相对误差 $e_r^{\ast }=\frac{x^{\ast }}{x}=\frac{x^{\ast }-x}{x}$
  • 在实际计算中 准确值 总是未知，通常取
    $$e_r^{\ast }=\frac{x^{\ast }}{x^{\ast }}=\frac{x^{\ast }-x}{x^{\ast }}$$
  • 假设误差不大时，两公式 取得 相对误差 相差极小，故忽略不计。
• 相对误差限: 相对误差的绝对值 上界
  • ${\epsilon }_r^{\ast }=\frac{{\epsilon }^{\ast }}{|x^{\ast }|}\ge |\frac{x^{\ast }-x}{x^{\ast }}|$

def 相对误差(x_star, x_true):
    return (x_star-x_true)/x_star

1.2. 有效数字

• 定义：$x^{\ast }$ 有 $n$ 位 有效数字。其表示为：
  $$\begin{gathered} x^{\ast }=\pm (a_0+a_1\times 10^{-1}+\cdots +a_{n-1}\times 10^{-(n-1)})\times 10^m \\ {a}_i(0\le i<n)\in I_{[0,9]},a_0\ne 0 \\ i,m\in I \end{gathered}$$
  $$|x-x^{\ast }|\le \frac{1}{2}\times 10^{m-n+1}$$

def 有效数字(star_x:str, x:float):
    string_x = str(star_x)
    if bool(re.match(r"\s*[+-]?([\d]+(\.[\d]*)?|\.[\d]+)(e[+-]?[\d]+)? *$", string_x)):
        float_x = float(star_x)
        if "e" not in string_x:
            string_x = "%e"%float_x
        M = int(string_x.split("e")[-1])
        error = np.abs(float_x - x) # |x - x^*|
        #print("误差为 %e" %误差(float_x, x))
        N = 0
        while (error <= (0.5 * 10**(M-N))):
            N += 1
        print(f"近似值 %.{N-1}e 有 {N} 位 有效数 -- 精确值：{x}。"%float_x)
        return M, N
    else:
        print("输入非数字！")

• 例：

有效数字("3.14", np.pi)
有效数字("3.1416", np.pi)
有效数字("8", 8.000033)
有效数字(2.7183, 2.7182818)
有效数字(-187.93, -187.93253)
M, N = 有效数字(0.037856, 0.03785551)

近似值 3.14e+00 有 3 位 有效数 -- 精确值：3.141592653589793。
近似值 3.1416e+00 有 5 位 有效数 -- 精确值：3.141592653589793。
近似值 8.0000e+00 有 5 位 有效数 -- 精确值：8.000033。
近似值 2.7183e+00 有 5 位 有效数 -- 精确值：2.7182818。
近似值 -1.8793e+02 有 5 位 有效数 -- 精确值：-187.93253。
近似值 3.7856e-02 有 5 位 有效数 -- 精确值：0.03785551。

• 相对误差限
  • 由 $x^{\ast }$ 定义, 可知
    $$\begin{array}{rl}{a}_0& \le |x^{\ast }|\times 10^{-m}<a_0+1\\ {\epsilon }_r^{\ast }& =\frac{{\epsilon }^{\ast }}{|x^{\ast }|}\le \frac{1}{2a_0}\times 10^{(m-n+1)-m}\\ & \le \frac{1}{2a_0}\times 10^{1-n}\end{array}$$
    当$x^{\ast }$ 具有 $n$ 位 有效数字时。
  • 反之 若 $x^{\ast }$ 相对误差限 ${\epsilon }_r^{\ast }\le \frac{1}{2(a_0+1)}\times 10^{1-n}$, 则 $x^{\ast }$ 至少具有 $n$ 位有效数字。
  • Q: 使 $\sqrt{1237}$ 的近似值相对误差小于 $0.1\%$，那么至少要取几位有效数字?
  • A: $\sqrt{1237}>\sqrt{1225}=35\to a_0=3$
    ${\epsilon }^{\ast }\le 0.166...\times 10^{-3}<10^{-3}=0.1\%$
    答：至少要取4位有效数字。
    验证：

有效数字('35.17',np.sqrt(1237))
print(f'相对误差为 {np.abs(相对误差(35.17,np.sqrt(1237)))} < 1.0e-3')

近似值 3.517e+01 有 4 位 有效数 -- 精确值：35.17101079013795。
相对误差为 2.8740123342252667e-05 < 1.0e-3

1.3. 误差估计

• 对两个近似数 $x_1^{\ast }$ 与 $x_2^{\ast }$ 的误差限分别为 $\epsilon (x_1^{\ast }),\epsilon (x_2^{\ast })$ ， 对于加减乘除运算的误差限满足不等式：

$$\begin{array}{rl}\epsilon (x_1^{\ast }\pm x_2^{\ast })& \le \epsilon (x_1^{\ast })+\epsilon (x_2^{\ast })\\ \epsilon (x_1^{\ast }\cdot x_2^{\ast })& \le |x_2^{\ast }|\cdot \epsilon (x_1^{\ast })+|x_1^{\ast }|\cdot \epsilon (x_2^{\ast })\\ \epsilon \left(\frac{x_1^{\ast }}{x_2^{\ast }}\right)& \le \frac{|x_2^{\ast }|\cdot \epsilon (x_1^{\ast })+|x_1^{\ast }|\cdot \epsilon (x_2^{\ast })}{|x_2^{\ast }{|}^2}，x_2^{\ast }\ne 0\end{array}$$

• 对于函数的误差计算，一般用泰勒展开进行估计，若 $f(x)$ 为可微方程：
  $$\begin{array}{rl}f(x)-f(x^{\ast })& =f^{\prime }(x^{\ast })(x-x^{\ast })+f^{\prime \prime }(\xi )(x-x^{\ast }{)}^2/2,\xi \in [x,x^{\ast }]\\ f(x)-f(x^{\ast })& =f^{\prime }(x^{\ast })e^{\ast }+f^{\prime \prime }(\xi ){e^{\ast }}^2/2\\ |f(x)-f(x^{\ast })|& \le |f^{\prime }(x^{\ast })|{\epsilon }^{\ast }+|f^{\prime \prime }(\xi )|{{\epsilon }^{\ast }}^2/2\end{array}$$
  • 一般 $f^{\prime \prime }(x^{\ast })$ 与 $f^{\prime }(x^{\ast })$ 比值相差不大，可忽略 $\epsilon$ 高阶项，计算函数的误差限 为：
    $$\epsilon (f(x^{\ast }))\approx |f^{\prime }(x^{\ast })|\epsilon (x^{\ast })$$
  • 拓展到多元函数 $A=f(x_1,x_2,\cdots ,x_n),A^{\ast }=f(x_1^{\ast },x_2^{\ast },\cdots ,x_n^{\ast })$
    $$\begin{array}{rl}e(A^{\ast })& =A^{\ast }-A\\ & \approx \sum _{k=1}^n(\frac{\partial f}{\partial x_k}{)}^{\ast }(x_k^{\ast }-x_k)\\ & \approx \sum _{k=1}^n\left(\frac{\partial A^{\ast }}{\partial x_k}\right)e_k^{\ast }\\ \epsilon (A^{\ast })& \approx \sum _{k=1}^n\left|\frac{\partial A^{\ast }}{\partial x_k}\right|\epsilon (x_k^{\ast })\\ {\epsilon }_r(A^{\ast })& =\frac{\epsilon (A^{\ast })}{|A^{\ast }|}\end{array}$$
• Q: 已测量 某操场 长宽为 $l^{\ast }=110m,w^{\ast }=80m$, 已知它们的误差限为 $|l-l^{\ast }|\le 0.2m,|w-w^{\ast }|\le 0.1m$，求面积的绝对误差限与相对误差限：
• A:
  $$\begin{array}{rl}s& =l\cdot w\\ {\frac{\partial s}{\partial l}}^{\ast }& =w^{\ast }=80m\\ {\frac{\partial s}{\partial w}}^{\ast }& =l^{\ast }=110m\\ \epsilon (s^{\ast })& \approx \left|{\frac{\partial s}{\partial l}}^{\ast }\right|\epsilon (w^{\ast })+\left|{\frac{\partial s}{\partial w}}^{\ast }\right|\epsilon (l^{\ast })\\ & \approx 80\times 0.2+110\times 0.1=27m^2\\ {\epsilon }_r(s^{\ast })& =\frac{\epsilon (s^{\ast })}{|s^{\ast }|}\approx \frac{27}{80\times 110}=0.31\%\end{array}$$
• 答：操场面积的绝对误差限为 $27m^2$，相对误差限为 $0.31\%$ 。