1、1第一章 绪 论本章以误差为主线，介绍了计算方法课程的特点，并概略描述了与算法相关的基本概念，如收敛性、稳定性，其次给出了误差的度量方法以及误差的传播规律，最后，结合数值实验指出了算法设计时应注意的问题1.1 引 言计算方法以科学与工程等领域所建立的数学模型为求解对象，目的是在有限的时间段内利用有限的计算工具计算出模型的有效解答。由于科学与工程问题的多样性和复杂性，所建立的数学模型也是各种各样的、复杂的. 复杂性表现在如下几个方面：求解系统的规模很大，多种因素之间的非线性耦合，海量的数据处理等等，这样就使得在其它课程中学到的分析求解方法因计算量庞大而不能得到计算结果，且更多的复杂数学模型没有分
2、析求解方法 这门课程则是针对从各种各样的数学模型中抽象出或转化出的典型问题，介绍有效的串行求解算法，它们包括(1) 非线性方程的近似求解方法；(2) 线性代数方程组的求解方法；(3) 函数的插值近似和数据的拟合近似；(4) 积分和微分的近似计算方法；(5) 常微分方程初值问题的数值解法；(6) 优化问题的近似解法；等等从如上内容可以看出，计算方法的显著特点之一是“近似” 之所以要进行近似计算，这与我们使用的工具、追求的目标、以及参与计算的数据来源等因素有关 计算机只能处理有限数据，只能区分、存储有限信息，而实数包含有无穷多个数据，这样，当把原始数据、中间数据、以及最终计算结果用机器数表示时就不
3、可避免的引入了误差，称之为舍入误差我们需要在有限的时间段内得到运算结果，就需要将无穷的计算过程截断，从而产生截断误差 如 的计算是无穷过程，当用!21e作为 的近似时，则需要进行有限过程的计算，但产生了!1!21nen截断误差 en2当用计算机计算 时，因为舍入误差的存在，我们也只能得到 的近似值ne ne，也就是说最终用 近似 ，该近似值既包含有舍入误差，也包含有截断误*e*差当参与计算的原始数据是从仪器中观测得来时，也不可避免得有观测误差 由于这些误差的大量存在，我们得到的只能是近似结果，进而对这些结果的“可靠性”进行分析就是必须的，它成为计算方法的第二个显著特点 可靠性分析包括原问题的适
4、定性和算法的收敛性、稳定性所谓适定性问题是指解存在、惟一，且解对原始数据具有连续依赖性的问题 对于非适定问题的求解，通常需要作特殊的预处理，然后才能做数值计算 在这里，如无特殊说明，都是对适定的问题进行求解对于给定的算法，若有限步内得不到精确解，则需研究其收敛性 收敛性是研究当允许计算时间越来越长时，是否能够得到越来越可靠的结果，也就是研究截断误差是否能够趋于零 对于给定的算法，稳定性分析是指随着计算过程的逐步向前推进，研究观测误差、舍入误差对计算结果的影响是否很大 对于同一类模型问题的求解算法可能不止一种，常希望从中选出高效可靠的求解算法 如我国南宋时期著名的数学家秦九韶就提出求 n 次多项
5、式值的如下快速算法011axxann；s； kntts),21(nk它通过 n 次乘法和 n 次加法就计算出了任意 n 次多项式的值 再如幂函数 可64x以通过如下快速算法计算出其值；xs；循环 6 次如上算法仅用了 6 次乘法运算，就得到运算结果算法最终需要在计算机上运行相应程序，才能得到结果，这样就要关注算法的时间复杂度(计算机运行程序所需时间的度量)、空间复杂度(程序、数据对存储空间需求的度量)和逻辑复杂度(关联程序的开发周期、可维护性以及可扩展性) 事实上，每一种算法都有自己的局限性和优点，仅仅理论分析是很不够的，大量的实际计算也非常重要，结合理论分析以及相当的数值算例结果才有可能选择
6、出适合自己关心问题的有效求解算法 也正因如此，只有理论分析结合实际计算才能真正把握准算法 31.2 误差的度量与传播一、误差的度量误差的度量方式有绝对误差、相对误差和有效数字定义 1.1 用 作为量 的近似，则称 为近似值 的绝对误差*x )(:*xe*x由于量 x 的真值通常未知，所以绝对误差不能依据定义求得，但根据测量工具或计算情况，可以估计出绝对误差绝对值的一个较小上界 ，即有(1.1)xe*)(称正数 为近似值 的绝对误差限，简称误差 这样得到不等式*x*工程中常用 x表示近似值 的精度或真值 x 所在的范围*x误差是有量纲的，所以仅误差数值的大小不足以刻划近似的准确程度 如量 (1.
7、2)mmcs 5012305.2315.0123为此，我们需要引入相对误差定义 1.2 用 作为量 的近似，称 为近似值 的相对误差 *xx)(:*xexr*x当 是 x 的较好近似时，也可以用如下公式计算相对误差*(1.3)*)(xer显然，相对误差是一个无量纲量，它不随使用单位变化 如式(1.2)中的量 s的近似，无论使用何种单位，它的相对误差都是同一个值 同样地，因为量 x 的真值未知，我们需要引入近似值 的相对误差限 ，*x)(*xr它是相对误差绝对值的较小上界 结合式(1.1)和(1.3)， 相对误差限可通过绝对误差限除以近似值的绝对值得到，即(1.4)*)(xr为给出近似数的一种表
8、示法，使之既能表示其大小，又能体现其精确程度，需引入有效数字以及有效数的概念定义 1.3 设量 的近似值 有如下标准形式x*xpnma21*.0(1.5) pmpnmaa 10101 其中 且 ，m 为近似值的量级 如果使不等式9,pia14(1.6)nmx102*成立的最大整数为 n，则称近似值 具有 n 位有效数字，它们分别是 、1a、 和 特别地，如果有 ，即最后一位数字也是有效数字，则称2ap是有效数*x从定义可以看出，近似数是有效数的充分必要条件是末位数字所在位置的单位一半是绝对误差限 利用该定义也可以证明，对真值进行“四舍五入”得到的是有效数 对于有效数，有效数字的位数等于从第一位
9、非零数字开始算起，该近似数具有的位数 注意，不能给有效数的末位之后随意添加零，否则就改变了它的精度例 1.1 设量 ，其近似值 ， ， 试回答这x14.3*x142.3*2x7*3x三个近似值分别有几位有效数字，它们是有效数吗？解 这三个近似值的量级 ，因为有m312*1 005.9.0x432 147852.3*31200.1x所以 和 都有 3 位有效数字，但不是有效数 具有 4 位有效数字，是有效*1 *2x数二、误差的传播这里仅介绍初值误差传播，即假设自变量带有误差，函数值的计算不引入新的误差 对于函数 有近似值 ，利用在点),(21nxfy ),(*2*1nxfy处的泰勒公式(Tay
10、lor Formula)，可以得到),(*2*1nx )(,()( *1*2* iii nxfye (1.7),(12ini nexf其中 ， 是 的近似值， 是 的绝对误差 式(1.7)iixf:*ii *ii ),21(ni表明函数值的绝对误差近似等于自变量绝对误差的线性组合，组合系数为相应的偏导数值 从式(1.7) 也可以推得如下函数值的相对误差传播近似计算公式5(1.8)(),()( *1*2* irini nr xeyxfye对于一元函数 ，从式(1.7)和(1.8) 可得到如下初值误差传播近似计算xf公式(1.9)()(*xefye(1.10)rr式(1.9)表明，当导数值的绝对值
11、很大时，即使自变量的绝对误差比较小，函数值的绝对误差也可能很大例 1.2 试建立函数 的绝对误差(限)、相对nnxxxfy 2121),(误差的近似传播公式，以及 时的相对误差限传播公式ii*0解 由公式(1.7) 和(1.8)可分别推得和的绝对误差、相对误差传播公式如下(1.11)niiini nxexfye 1*1*2* )()(,()( (1.12)iiririir yyx*),进而有niiniinii xeye1*1*1* )()()()( 于是有和的绝对误差限近似传播公式niix1*)()(当 时，由式(1.3)推得相对误差限的近似传播公式niix1*0 )(max)(ax )()(
12、 *1*1 1*1*1* irninirni niirniniiriniir yyxy 例 1.3 使用足够长且最小刻度为 1mm 的尺子，量得某桌面长的近似值mm，宽的近似值 mm (数据的最后一位均为估计值) 试求3.104*a 8.704*b桌子面积近似值的绝对误差限和相对误差限 解 长和宽的近似值的最后一位都是估计位，尺子的最小刻度是毫米，故有误差限mm， mm5.0)(*a5.0)(*b6面积 ，由式(1.7) 得到近似值 的绝对误差近似为abS*baS)()()(*ee进而有绝对误差限mm25.104.3145.0874)()(* a相对误差限%3.10* Sr1.3 数值实验与算
13、法性能比较本节通过几个简单算例说明解决同一个问题可以有不同的算法，但算法的性能并不完全相同，他们各自有自己的适用范围，并进而指出算法设计时应该注意的事项算例 1.1 表达式 ，在计算过程中保留 7 位有效数字，研究)1(1xx对不同的 x，两种计算公式的计算精度的差异说明 1：Matlab 软件采用 IEEE 规定的双精度浮点系统，即 64 位浮点系统，其中尾数占 52 位，阶码占 10 位，尾数以及阶码的符号各占 1 位 机器数的相对误差限(机器精度)eps=2 52 2.22044610 16 ，能够表示的数的绝对值在区间(2.225073910308 ，1.79769310 308)内，
14、该区间内的数能够近似表达，但有舍入误差，能够保留至少 15 位有效数字 其原理可参阅参考文献2, 4分析算法 1： 和算法 2： 的误差时，精确解用1)(xy )1()xy双精度的计算结果代替 我们选取点集 中的点作为 x，比较两种方法误差301i的差异从图 1.1 可以看出，当 x 不是很大时，两种算法的精度相当，但当 x 很大时算法 2 的精度明显高于算法 1 这是因为，当 x 很大时， 和 是相近数，用x1算法 1 进行计算时出现相近数相减，相同的有效数字相减后变成零，于是有效数字位数急剧减少，自然相对误差增大 这一事实也可以从误差传播公式(1.12)分析出 鉴于此，算法设计时，应该避免
15、相近数相减在图 1.2 中我们给出了当 x 接近 时，两种算法的精度比较，其中变量 x 依次1取为 从图中可以看出两种方法的相对误差基本上都为 ，因而二301ii 710者的精度相当7图 1.1 算例 1.1 中两种算法的相对误差图（ ）x图 1.2 算例 1.1 中两种算法的精度比较 )1(x算例 1.2 试用不同位数的浮点数系统求解如下线性方程组 2310.1x说明 2：浮点数系统中的加减法在运算时，首先按较大的阶对齐，其次对尾数实施相应的加减法运算，最后规范化存入计算机 算法 1 首先用第一个方程乘以适当的系数加至第二个方程，使得第二个方程的 的系数为零，这时可解出 ；其次将 带入第一个
16、方程，进而求得 (在第x2x2x1x三章中称该方法为高斯消元法) 当用 4 位和 7 位尾数的浮点运算实现该算法，分别记之为算法 1a 和算法 1b8算法 2 首先交换两个方程的位置，其次按算法 1 计算未知数 (第三章中称其为选主元的高斯消元法) 当用 4 位和 7 位尾数的浮点运算实现该算法，分别记之为算法 2a 和算法 2b方程组的精确解为 ， ，用不同的算法计算.25018.1x .4987.02x出的结果见表 1.1 表 1.1 对算例 1.2 用不同算法的计算结果比较算例1.2*1x)(*1xr*2x)(*2xr算法 1a 0.0000 0.10101 0.5000 0.25107
17、算法 2a 0.2500 0.75107 0.5000 0.25107算法 1b 0.2600000 0.40101 0.4999987 0.10106算法2b0.2500020 0.50108 0.5000000 0.25107对于算例 1.2，表中的数据表明，当用 4 位尾数计算时，算法 1 给出错误的结果，算法 2 则给出解很好的近似 这是因为在实现算法 1 时，需要给第一个方程乘以 加至第二个方程，从而削去第二个方程中 的系数，但在计01./ x算 的系数时需做如下运算2x(1.13)6616 03.04.3.04.3. 对上式用 4 位尾数进行计算，其结果为 因为舍入误差，给相对较大
18、的数加以相对较小的数时，出现大数“吃掉”小数的现象 计算右端项时，需做如下运算 (1.14)6616 102.0202.10.2 同样出现了大数吃小数现象，其结果为 这样，得到的变形方程组6261114.xx中没有原方程组中第二个方程的信息，因而其解远偏离于原方程组的解 该算法中之所以出现较大数的原因是因为运算 ，因而算法设计中尽可能0./避免用绝对值较大的数除以绝对值较小的数 其实当分子的量级远远大于分母的量级时，除法运算还会导致溢出，计算机终止运行虽从单纯的一步计算来看，大数吃掉小数，只是精度有所损失，但多次的大数吃小数，累计起来可能带来巨大的误差，甚至导致错误 例如在算法 1a 中出9现
19、了两次大数吃小数现象，带来严重的后果 因而尽可能避免大数吃小数的出现在算法设计中也是非常必要的当用较多的尾数位数进行计算，舍入误差减小，算法 1 和 2 的结果都有所改善，算法 1 的改进幅度更大些算例 1.3 计算积分 有递推公式 ，已知105dxIn ),(51nIIn 采用 IEEE 双精度浮点数，分别用如下两种算法计算 的近似值56ln0I 30算法 1 取 的近似值为 ，按递推公式 计0I 6793.1825*0I *1*5nnII算 *30I算法 2 因为 ，取 的近似)19(5)39(610391039 dxIdx 39I值为，按递推公式 计算.458014*39I *1nnII
20、*30I算法 1 和算法 2 的计算结果见表 1.2 误差绝对值的对数图见图 1.3表 1.2 算例 1.3 的计算结果算法 1 算法 2n *InI*n *InI*1 8.8392e-002 1.9429e-016 39 4.5833e-003 3.9959e-0042 5.8039e-002 9.8532e-016 38 4.2115e-003 7.9919e-0053 4.3139e-002 4.9197e-015 37 4.4209e-003 1.5984e-0054 3.4306e-002 2.4605e-014 36 4.5212e-003 3.1967e-0065 2.8468e
21、-002 1.2304e-013 35 4.6513e-003 6.3935e-0076 2.4325e-002 6.1520e-013 34 4.7840e-003 1.2787e-007 33 4.9255e-003 2.5574e-00825 1.1740e+001 1.1734e+001 32 5.0755e-003 5.1148e-00926 -5.8664e+001 5.8670e+001 31 5.2349e-003 1.0230e-00927 2.9336e+002 2.9335e+002 30 5.4046e-003 2.0459e-01028 -1.4667e+003 1.
22、4668e+003 29 7.3338e+003 7.3338e+003 30 -3.6669e+004 3.6669e+00410图 1.3 算例 1.3 用不同算法计算结果的误差绝对值的对数图从表 1.2 中的计算结果可以看出，算法 1 随着计算过程的推进，绝对误差几乎不断地以 5 的倍数增长，即有 0*2*21* 555 IIII nnnn 成立 对于逐步向前推进的算法，若随着过程的进行，相对误差在不断增长，导致产生不可靠的结果，这种算法称之为数值不稳定的算法 对于算法 1 绝对误差按 5 的幂次增长，但真值的绝对值却在不断变小且小于 1，相对误差增长的速度快于 5 的幂次，导致产生错误
23、的结果，因而算法 1 数值不稳定，不能使用 而算法 2 随着计算过程的推进，绝对误差几乎不断地缩小为上一步的 1/5，即有mnmnnn IIII 5/5/5/ *2*21* 成立 绝对误差不断变小，真值的绝对值随着过程向前推进却在变大，这样相对误差也越来越小，这样的方法称之为数值稳定的算法 算法 1 和算法 2 的误差对数示意图见图 1.3 这个算例告诉我们应该选用数值稳定的算法知识结构图11算 法 设 计 要 点数 值 方 法 的 稳 定 性数 值 方 法 的 收 敛 性算 法 多 元 函 数一 元 函 数传 播 有 效 数 字相 对 误 差 （ 限 ）绝 对 误 差 （ 限 ）度 量 截
24、断 误 差舍 入 误 差误 差 的 产 生误 差误 差 与 算 法习题一1 已知有效数 ， ， 试给出各个近似值的105.3*x4*2105.x10.*3x绝对误差限和相对误差限，并指出它们各有几位有效数字2 证明当近似值 是 x 的较好近似时，计算相对误差的计算公式 和* x*相差一个和 同阶的无穷小量*x23 设 x 的近似值 具有如式(1.5)的表示形式，试证明*x1) 若 具有 n 位有效数字，则相对误差 ；nraxe1*02)(2) 若相对误差 ，则 至少具有 n 位有效数字nrae11*0)(2)4 试建立二元算术运算的绝对误差限传播近似计算公式5 试建立如下表达式的相对误差限近似
25、传播公式，并针对第 1 题中数据，求下列各近似值的相对误差限1) ； 2) ； 3) *321*xy3*2xy*32*3/xy6 若例题 1.3 中使用的尺子长度是 80mm，最小刻度为 1mm，量得某桌面长的近似值 mm，宽的近似值 mm 试估计桌子长度、宽度的绝.04a 8.704b对误差限，并求用该近似数据计算出的桌子面积的绝对误差限和相对误差限7 改变如下计算公式，使其计算结果更为精确1) 且,cos1x12) 1,ln)ln()(ln NNdN123) 1,133xx8 (数值试验 )试通过分析和数值试验两种手段，比较如下三种计算 近似值算法1e的可靠性算法 1 ；mne0!)(算法 2 ；1算法 3 ；101)!(mne9 (数值试验 )设某应用问题归结为如下递推计算公式， ，72.80y251ny,1在计算时 取为具有 5 位有效数字的有效数 试分析近似计算公式*c的绝对误差传播以及相对误差传播情况，并通过数值实验验证 *1*cyn(准确值可以用 IEEE 双精度浮点运算结果代替)，该算法可靠可用吗？
