《数值计算方法》系列总目录
第一章 误差序列实验
第二章 非线性方程f(x)=0求根的数值方法
第三章 CAD模型旋转和AX=B的数值方法
第四章 插值与多项式逼近的数值计算方法
第五章 曲线拟合的数值方法
第六章 数值微分计算方法
第七章 数值积分计算方法
第八章 数值优化方法
数值计算往往关注各种数值计算方法的误差、稳定性、收敛性、计算工作量、存贮量和自适应性,因为这些因素往往关乎于进行科学计算最终得出的结果的正确性以及效率。
因此,基于计算机的发展及其在各技术科学领域的应用推广与深化,新的计算性学科分支纷纷兴起,数值计算方法广泛应用于计算力学、计算物理、计算电磁学、高性能并行计算、数字图像处理等一系列工程学科中。目前常用的一些数值计算软件,例如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=34pn−1−31pn−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=310qn−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++)