误差序列实验——《数值计算方法》

《数值计算方法》系列总目录

第一章 误差序列实验
第二章 非线性方程f(x)=0求根的数值方法
第三章 CAD模型旋转和AX=B的数值方法
第四章 插值与多项式逼近的数值计算方法
第五章 曲线拟合的数值方法
第六章 数值微分计算方法
第七章 数值积分计算方法
第八章 数值优化方法


第一章


系列前言

数值计算方法在计算机科学中的应用

   数值计算方法属于计算数学的一个分支,而计算数学又是属于数学科学的一个分支学科。计算数学是现代数学的一个分支,基于计算机技术的迅猛发展而兴起的一个学科。计算数学主要研究内容为数值计算方法,即研究能够使用计算机求解实际问题的算法与数学原理。通俗来讲,就是为计算机“量身打造”求解实际问题的数学公式定理,其中最明显的特点即为迭代算法与差分算法。
   计算数学几乎与数学科学的一切分支有联系。一方面,它的根本原理来源于数学领域的公式定理,并在其基础上发展了可以计算机实现的高性能算法,使得计算机技术可以广泛的运用于现代科学研究,省去了大量的繁琐的人工复杂计算。另一方面,绝大多数现代科学于工程研究都需要大量的计算,比如CAD技术、电磁仿真技术、材料性质的数理方程计算等等。可见计算数学是架起计算机与实际问题需求的一个重要的桥梁,而优越的数值计算方法就显得尤为重要。因此,计算数学的内容十分广泛。
   面向现实问题来看,目前全世界的各个科技行业与研究领域都急需数值算法方面的人才,这与计算数学与计算机技术的强大的优越性密不可分。在实际工程开发中,都需要大量的科学计算,然而由于计算机的性能有限以及实际问题的刚性需求,必须对算法进行改进以求得更加优越、复杂度更小、通用性更高的算法,因此需要开发出更加强大的数值算法。其中经典的可实际工程化的算法包括FDTD算法、最速下降法等。

数值计算往往关注各种数值计算方法的误差、稳定性、收敛性、计算工作量、存贮量和自适应性,因为这些因素往往关乎于进行科学计算最终得出的结果的正确性以及效率。

  因此,基于计算机的发展及其在各技术科学领域的应用推广与深化,新的计算性学科分支纷纷兴起,数值计算方法广泛应用于计算力学、计算物理、计算电磁学、高性能并行计算、数字图像处理等一系列工程学科中。目前常用的一些数值计算软件,例如MATLAB、Mathematica、CST、Comsol和HFS等都是基于各个专业领域的高性能数值方法开发的,并在现代科技行业内大放异彩。
  综上所述,数值计算方法本质来源于计算数学的蓬勃发展,并与其始终进行频繁的交互,在此前提下,再基于现代计算机科学技术,使得数值计算方法成为联系数学与计算机两大学科的重要枢纽,并应用于包括计算数学和计算机科学在内的一大批工程计算与开发问题,成为促进现代工程与科学飞速发展的不可或缺的一大部分。


一、算法原理

   由于计算机进行数值运算的精度限制以及算法的设计,大量次的运算会导致运算结果与期望的真实值之间产生误差。计算机进行运算时的误差主要来源于舍入误差,而人为设计的算法的误差主要来源于截断误差。

  在计算机进行连续计算过程中,会产生误差传播。考虑p和q的加法运算,他们的近似值分别是和,误差分别是和。从p= p ^ \hat p p^​+ ε p {\varepsilon _p} εp​和p= p ^ \hat p p^​+ ε p {\varepsilon _p} εp​开始,它们的和为
    p + q = ( q ^ + ε q ) + ( p ^ + ε p ) = ( q ^ + p ^ ) + ( ε q + ε p ) p + q = (\hat q + {\varepsilon _q}) + (\hat p + {\varepsilon _p}) = (\hat q + \hat p) + ({\varepsilon _q} + {\varepsilon _p}) p+q=(q^​+εq​)+(p^​+εp​)=(q^​+p^​)+(εq​+εp​)
  因此对于加法运算,和的误差是每个加数的误差的和。
  在乘积计算过程中,误差的传播更为复杂。其中,乘积的表达式如下:
p q = ( p ^ + ε p ) ( q ^ + ε q ) = p ^ q ^ + p ^ ε q + q ^ ε p + ε q ε p pq = (\hat p + {\varepsilon _p})(\hat q + {\varepsilon _q}) = \hat p\hat q + \hat p{\varepsilon _q} + \hat q{\varepsilon _p} + {\varepsilon _q}{\varepsilon _p} pq=(p^​+εp​)(q^​+εq​)=p^​q^​+p^​εq​+q^​εp​+εq​εp​
  由公式可知,pq的误差与和的取值情况相关,因此可以考虑pq乘积的相对误差。有相对误差的公式可以计算得:
R p q = p q − p ^ q ^ p q = p ^ ε q p q + q ^ ε p p q + ε q ε p p q ≈ ε q q + ε p p + 0 = R q + R p {R_{pq}} = \frac{{pq - \hat p\hat q}}{{pq}} = \frac{{\hat p{\varepsilon _q}}}{{pq}} + \frac{{\hat q{\varepsilon _p}}}{{pq}} + \frac{{{\varepsilon _q}{\varepsilon _p}}}{{pq}} \approx \frac{{{\varepsilon _q}}}{q} + \frac{{{\varepsilon _p}}}{p} + 0 = {R_q} + {R_p} Rpq​=pqpq−p^​q^​​=pqp^​εq​​+pqq^​εp​​+pqεq​εp​​≈qεq​​+pεp​​+0=Rq​+Rp​
  发现乘积pq的相对误差大致等于和的相对误差之和。
  因此如果数据具有初始误差,那么会通过一系列的计算进行传播,因此必须设计算法尽量减少初始误差传播对最终结果产生的影响较小,这样的算法称为稳定算法,否则称为不稳定算法。

二、实验内容及核心算法代码

   为研究不同计算方法对误差传播以及最终结果的影响,设计可进行连续运算的算法并人为引入一定的初始误差并比较每次运算的结果与真实值的误差,以此来仿真模拟初始误差的传播特性以及误差的敛散性。

  首先设计用无限精度算法结合下列3个算法可递归生成序列的各项值。
r 0 = 1 , a n d r n = r n − 1 / 3 {r_0} = 1,{\rm{ and }}{r_n} = {r_{n - 1}}/3 r0​=1,andrn​=rn−1​/3
p 0 = 1 , p 1 = 1 3 , a n d p n = 4 3 p n − 1 − 1 3 p n − 2 {p_0} = 1,{p_1} = \frac{1}{3},{\rm{ and }}{p_n} = \frac{4}{3}{p_{n - 1}} - \frac{1}{3}{p_{n - 2}} p0​=1,p1​=31​,andpn​=34​pn−1​−31​pn−2​
q 0 = 1 , q 1 = 1 3 , a n d q n = 10 3 q n − 1 − q n − 2 {q_0} = 1,{q_1} = \frac{1}{3},{\rm{ and }}{q_n} = \frac{{10}}{3}{q_{n - 1}} - {q_{n - 2}} q0​=1,q1​=31​,andqn​=310​qn−1​−qn−2​
然后,通过设置 r 0 = 0.99996 {r_0} = 0.99996 r0​=0.99996, p 0 = 1 , p 1 = 0.33332 {p_0} = 1,{p_1} = 0.33332 p0​=1,p1​=0.33332, q 0 = 1 , q 1 = 0.33332 {q_0} = 1,{q_1} = 0.33332 q0​=1,q1​=0.33332用来表示序列的近似值,通过比较三个序列的每一项的值与真实值的误差以及其相对误差,来模拟在连续运算中初始误差的传播情况。
以下为C++核心算法代码:

#pragma once
#if 1
#include <iostream>
#include <cmath>
#include <iomanip>
using namespace std;
void ErrSequenceList(float v0, float v1)
{
	//定义每个算法结果序列
	const int NUM = 15;
	float x[NUM], r[NUM], p[NUM], q[NUM];
	x[0] = 1, x[1] = x[0] / 3.;
	r[0] = 0.99996,r[1]=r[0]/3;
	p[0] = q[0] = v0;
	p[1] = q[1] = v1;
	//计算环节
	for (int i = 2; i < NUM; i++)
	{
		x[i] = 1 / powf(3, i);
		r[i] = r[i - 1] / 3;
		p[i] = 4 * p[i - 1] / 3 - p[i - 2] / 3;
		q[i] = 10 * q[i - 1] / 3 - q[i - 2];
	}
	//格式化输出
	cout << "sequence of three recursive value:" << endl;
	cout << setw(20) << "n" << setw(20) << "xn";
	cout << setw(20) << "rn" << setw(20) << "pn" << setw(20) << "qn" << endl;
	cout << setiosflags(ios::fixed) << setprecision(10);
	for (int i = 0; i < NUM; i++)
	{
		cout << setw(20) << i << setw(20) << x[i];
		cout << setw(20) << r[i] << setw(20) << p[i] << setw(20) << q[i] << endl;
	}
	cout << "error sequence of three recursive value:" << endl;
	cout << setw(20) << "n" << setw(20) << "|x[n]-r[n]|";
	cout << setw(20) << "|x[n]-p[n]|" << setw(20) << "|x[n]-q[n]|" << endl;
	for (int n = 0; n < NUM; n++)
	{
		cout << setw(20) << n << setw(20) << fabs(x[n] - r[n]);
		cout << setw(20) << fabs(x[n] - p[n]) << setw(20) << fabs(x[n] - q[n]) << endl;
	}
}
int main()
{
	float v0 = 1, v1 = 1. / 3;
	ErrSequenceList(v0, v1);
	return 0;
}
#endif

三、计算结果及分析

误差序列实验——《数值计算方法》

                              图1 不稳定的误差序列{Xn-Qn}(python matplotlib绘图库)


四、结论

  在相同的初始误差情况下,不同的算法会导致在连续计算中误差的传播情况产生不同的结果,即收敛与发散。在进行数值计算时应当考虑算法对误差传播的影响,即相对误差的敛散性与稳定性,尽量采取稳定算法。


声明:本系列《数值计算方法》所有章节源代码均为作者编写,大家可以自行下载,仅供学习交流使用,欢迎大家评论留言或者私信交流,请勿用于任何商业行为。

各章算法源代码(C++)

上一篇:MATLAB卡尔曼滤波-实例


下一篇:线性回归详解