西安交通大学计算方法C讲义

1、计算方法（C）目 录第1章 绪论1.1 数值计算1.2 数值方法的分析1.2.1 计算机上数的运算1.2.2 算法分析第2章 线性代数方程组 2.1 Gauss消去法2.1.1 消去法2.1.2 主元消去法2.2 矩阵分解2.2.1 Gauss消去法的矩阵意义2.2.2 矩阵的LU分解及其应用2.2.3 其他类型矩阵的分解2.2.4 解三对角矩阵的追赶法2.3 线性方程组解的可靠性2.3.1 向量和矩阵范数2.3.2 残向量与误差的代数表征2.4 解线性方程组解的迭代法2.4.1 基本迭代法2.4.2 迭代法的矩阵表示2.4.3 收敛性第3章 数据近似3.1 多项式插值3.1.1 插值多项式3

2、.1.2 Lagrange插值多项式3.1.3 Newton插值多项式3.1.4 带导数条件的插值多项式3.1.5 插值公式的余项3. 2 最小二乘近似最小二乘问题的法方程 正交化算法第4章 数值微积分4.1 内插求积，Newton-Cotes公式4.1.1 Newton-Cotes公式4.1.2 复化求积公式4.1.3 步长的选取4.1.4 Romberg方法4.1.5 待定系数法4.2 数值微分4.2.1 插值公式方法4.2.2 Taylor公式方法 (待定系数法)4.2.3 外推法第5章 非线性方程求解5.1 解一元方程的迭代法5.1.1 简单迭代法5.1.2 Newton法5.1.3

3、割线法5.1.4 区间方法5.2 收敛性问题5.2.1 简单迭代不动点5.2.2 收敛性的改善5.2.3 Newton法的收敛性5.2.4 收敛速度第1章 绪论1.1数值计算现代科学的发展，已导致科学与技术的研究从定性前进到定量，尤其是现代数字计算机的出现及迅速发展，为复杂数学问题的定量研究与解决，提供了强有力的基础。通常我们面对的理论与技术问题，绝大多数都可以从其物理模型中抽象出数学模型，因此，求解这些数学模型已成为我们面临的重要任务。一、 本课程的任务：寻求解决各种数学问题的数值方法如何将高等数学的问题回归到初等数学（算术）的方法求解了解计算的基础方法，基本结构（否则只须知道数值软件）并研

4、究其性质。立足点：面向数学解决数学问题面向计算机利用计算机作为工具充分发挥计算机的功能，设计算法，解决数学问题 例如：迭代法、并行算法二、 问题的类型1、离散问题：例如，求解线性方程组 从离散数据：矩阵A和向量b，求解离散数据x；2、连续问题的离散化处理：例如，数值积分、数值微分、微分方程数值解；3、离散问题的连续化处理：例如，数据近似，统计分析计算；1.2数值方法的分析在本章中我们不具体讨论算法，首先讨论算法分析的基础误差。一般来讲，误差主要有两类、三种（对科学计算）：1）公式误差“截断误差”，数学计算，算法形成主观（人为）：数学问题-数值方法的转换，用离散公式近似连续的数学函数进行计算时，

5、一般都会发生误差，通常称之为“截断误差”； 以后讨论2）舍入误差及输出入误差计算机，算法执行客观（机器）：由于计算机的存储器、运算器的字长有限，在运算和存储中必然会发生最末若干位数字的舍入，形成舍入误差；在人机数据交换过程中，十进制数和二进制数的转换也会导致误差发生，这就是输入误差。这两种误差主要是由于计算机的字长有限，采用浮点数系所致。首先介绍浮点数系计算机上的运算浮点运算面向计算机设计的算法，则先要讨论在计算机上数的表示。科学记数法浮点数：约定尾数中小数点之前的数全为零，小数点后第一个数不能为零。目前，一般计算机都采用浮点数系，一个存储单元分成首数和尾数：×××

6、;××××首数 尾数(位)其中首数存放数的指数（或“阶”）部分，尾数存放有效数字。对于b进制，尾数字长为t位的浮点数系中的（浮点）数，可以用以下形式表示：此处，指数（称为阶）限制在范围内。以下记实数系中的实数为，它在浮点数系中对应的浮点数记为进制，尾数位数，阶的范围。几乎所有近代计算机都采用“二进制”（即）：位、字节和字分别是指位数不同的二进制数。例如十进制转换二进制10000 000120000 001040000 010080000 100090000 1001100000 101027位是一个二进制数（即0或1）；字节是8个二进制数字；上表的最后一

7、列是字节。单精度浮点数（single precision）按32位存储，双精度浮点数（double precision）按64位存储。精度用于指明每个浮点数保留多少位以及尾数和阶数各分配多少位。单精度浮点数的尾数为23位、阶数为8位；双精度浮点数的尾数为53位（包含符号位）、阶数为11位（包含符号位）。双精度浮点数的等价二进制数如下所示：浮点数的特点：1、 实数转换到浮点数浮点化，缺点：总会产生误差（除极个别的情况：）按 四舍五入，绝对误差：（举例），优点：浮点化产生的相对误差有界（与数字本身的数量级无关） 注：设实数，则按进制可表达为：按四舍五入的原则，当它进入浮点数系时，若，则若，则 对第

8、一种情况：对第二种情况：就是说总有： 另一方面，由浮点数要求的 ， 有，将此两者相除，便得2、 每一个浮点数系的数字有限： 3、 浮点数系中的运算非自封闭,（因为数字有限、尾数字长有限、指数数字有限等）例：在中，运算和,的结果显然已不在此浮点数系内，而或也不在此浮点数系内，需结果才在此浮点数系内。浮点运算应注意：1） 避免产生大结果的运算，尤其是避免小数作为除数参加运算；2） 避免“大”“小”数相加减；3） 避免相近数相减，防止大量有效数字损失；4）尽可能简化运算步骤，减少运算次数。原因：记，由，可得： （“。”表示任意一种四则运算）此处 是由机器字长（实质上是尾数字长的大小）确定的常数（它反

9、映了实际运算的精度）。显然，若要求运算的舍入误差小，应使运算结果（如：）较小。尤其是小分母运算：，小大误差。其次，当浮点数系中两个数量级相差较大的数相加（或减），注意到由于浮点数系中数字字长的有限性，可能导致“大数吃小数”。例如，中，则似乎没有参加运算。 第三，同样，由于浮点数系中数字字长的有限性，当两个相近数相减时：例如，在中，两数相减：，计算结果仅剩2位有效数字，而原来参加运算的数字有8位有效数字，这将严重影响最终计算结果的精度。 算法分析作为一个可用的算法，必须考虑其效率和可靠性，定义：计算机完成一个乘法和一个加法，称为一个浮点运算（记为flop）；注：由于计算机在运算时，加（减）法所耗

10、时间远少于乘（除）法，所以通常只须计算乘法的次数，因此也有“一个算法需要多少个乘法”这一提法。1、计算效率可计算性（计算复杂性空间、时间）例：解线性方程组的Grammar 方法：，其中是方程组系数矩阵对应的行列式，而则是以右端向量替代的第列所得矩阵对应的行列式。由线性代数知识可知，若，则，其中是由变换到所需的置换次数。可见每计算一个行列式，需要个浮点运算；因此，按Grammar方法解方程组约需 个浮点运算。当时，用一个运算速度为的计算机进行求解，约需年（日前报道我国计算机已达到秒，这仍需近10年）。而n=20的方程组应该说是一个小型的方程组。因此Grammar方法是一个不能接受的算法，它缺乏可

11、计算性。第二章将介绍的Gauss消去法和迭代法就有较高的效率，具有很好的可计算性。2、计算可靠性作为算法，除了考虑其效率外，必须重视可靠性，它包括两方面：问题的性态 和 方法的稳定性l 问题的性态所计算的问题当原始数据发生小扰动时，问题的解一般也发生扰动。问题的性态问题的解对原始数据发生变化的敏感性。原始数据小扰动问题解 例：线性方程组： 的解是：若将方程组的系数改写成具有2位有效数字的小数： 的解则变成：；这是一个典型的病态方程组。一般：由原始数据 计算结果， 扰动后的数据 计算结果，若对问题存在常数m，满足关系式：或 则称（相对误差之比的上界）m为该问题的条件数，记作 ；由微积分中值定理知

12、识容易计算出，当时， 。稍后我们在第二章将对线性方程组的性态作进一步的分析。l 算法的数值稳定性：例：计算；解：由微积分知识，可得计算公式， ，我们将准确值与计算结果列表如下：方法准确值.6321.3679.2642.2073.1709.1455.1268.1124.1010.6321.3679.2642.2074.1704.1480.1120.2160-.7280.6320.3680.2643.2073.1709.1455.1268.1124.1010由上表可见，方法中，原始步的误差，随着计算步数的增加被严重地放大，特别是竟变成负数（注意：被积函数是非负函数），而方法则相反；这是因为方法中，

13、若前步有误差：，则，说明的误差，经一步计算后被扩大了倍，随着的增大，误差将被大大地扩大；而通过同样的分析可知：方法中的误差则被缩小倍。算法的数值稳定性：算法对初始误差导致的最终误差的可控性，如果最终误差被有效地控制，则称算法是稳定的，否则就是不稳定的。第二章 线性方程组求解线性方程组：其中2.1消去法 消去法设计方法原则：面向计算机（事先未知元素，编制程序），例：基本思想：降维N维问题转化为N-1维问题逐次降维，依次进行消去过程对方程组（2.1）由上而下逐步消去对角元 以下的 ，使之转化为如下等价形式，达到目标：从而，可进入回代过程，再由下而上，逐步解得 ：这儿的主元对问题（2.1）设，目标：

14、将第2n方程的的系数转化为0；方法：“第个方程”-“第1个方程”（，得到现在只须关心由第2n方程形成的n-1维方程组，以后可仿此进行。消去：第步：设，以第行为基准，消去以下各行中列以下的，令 施行：第行第行新的第行，元素变化：（），经过步消去（注意：以下无元可消），得到式。注意：每计算1个仅需1个浮点运算回代：第一步 ， 第二步 ，第步 运算量：消元：n元方程组：第1步消元：对第行，共n-1行；每1行：计算，1个乘除法（或“浮点运算”），计算新的共 个乘除法（或“浮点运算”），第1次元共个乘除法，因此，消元的运算量； 回代：第1步，求需要1个乘除法（或“浮点运算”），即：第2步求用2个乘除法（

15、或“浮点运算”），一般，第步求用个浮点运算，因此，回代的运算量；消去法 的 总运算量： 。例2-1：解线性方程组解：利用增广矩阵（因为线性方程组的解仅与系数与右端常数项相关）例： ；解；解：；解 ；解；解；2.1.2 （列）主元消去法算法中，若第步：或则按照原来的简单消去法算法可能无法执行，也可能出现大误差：例：在浮点数系中计算；方程组 准确解：解：按消去法解：误差大，原因：若有误差 则，则演变为；这说明前步各元素的误差均被放大倍。克服方法，将按绝对值最大的元素交换到主元位置，使，使前步的误差不再被放大。消元过程中，第步消元以第行为基准，消去其后的个方程的（系数矩阵第列以下的元素），消去法要求

16、主元。为避免出现作为主元，并使前步的误差不被放大，应使，为此应使：，通常采用按列选主元的列主元消去法：若，便将第行与第行交换，使与交换位置，使新的第行执行在原始消去法中的角色，保证将作为除数的主元，从而，重复前述的消去过程。列主元消去法的效果：1. 算法稳定，即计算误差能被有效控制；2. 当矩阵的行列式时，算法总可以执行完成；注：若矩阵是对称正定或严格对角占优，则不选主元，消去法也是稳定的；矩阵严格对角占优的定义：对角元的绝对值大于该行其他元的绝对值之和，即；问题：为什么系数矩阵严格对角占优时，不选主元的消去法也是稳定的？为保证消去法是稳定的，计算应如何进行？注意：第一步消元 ，由于占优：，它

17、的效果与主元消去法的一样，误差不会被放大。但因此算法也应适当改变，应先对第一行计算此值，然后从第2列起逐列从上到下计算。且第一步消元后生成的右下方阶矩阵仍是严格对角占优，以第二列为例：又即新的第二行的对角元绝对占优，其他也可同样证明。例2-2：列主元消去法求解例2-1同样得到原方程组的解 ；2.2 矩阵分解 Gauss消去法的矩阵意义矩阵的三角分解若将消去法解方程过程中产生的作为一个单位下三角矩阵其对角元为1，对角线上的元素均为0的对应元素；将消去法消元过程最后得到的上三角矩阵对角线以下元素均为0记作或（仅考虑方程的系数矩阵部分）；如例2.1，再将与或相乘：或显然相乘得到的正是原始所求解的方程

18、的增广矩阵，而相乘得到的正是原始所求解的方程的系数矩阵。反过来，一个矩阵也可以通过求解的过程获得这样的“分解”，其中为单位下三角矩阵，为上三角矩阵。这样将一个矩阵分解成一个下三角矩阵和一个上三角矩阵的乘积形式，称为矩阵的三角（Doolittle）分解。1） 消去法的矩阵意义： 以例2-1说明与以下矩阵乘法是一样的：可见，一般的第一步消元是：若记第步消元后形成的矩阵为，特别：则第步消元是将以下初等矩阵左乘矩阵形成因此，Gauss消去过程从矩阵运算的角度来看，是其中 为上三角矩阵，为向量：矩阵分解注意到初等矩阵具有性质：及因此，我们有 根据矩阵乘积的性质，有 这就是基本的矩阵三角分解Doolitt

19、le分解将矩阵分解为单位下三角矩阵与上三角矩阵的乘积。从矩阵乘法与行列式的关系可知，若矩阵A三角分解得到A=LU，则有：det A = det L * det U，由于下三角矩阵或上三角矩阵的行列式是全部对角元的乘积，因此可利用三角分解求矩阵对应的行列式的值。例如：上述例2.1方程组系数矩阵对应的行列式的值是-1，下例2.3 对称正定矩阵对应的行列式的值为144。当系数矩阵为的方程组可以顺利求解（不必选主元）时，解的过程显然是唯一的，因此它的Doolittle分解必然也是唯一的，从而可以通过待定系数的方法获得该矩阵的Doolittle分解（通常称为“直接分解”或“紧凑格式”）。 其他三角分解除

20、了上述的单位下三角矩阵与上三角矩阵的乘积形式以外，还有其他类型的分解，例如下三角矩阵和单位上三角矩阵的乘积，单位下三角矩阵与对角矩阵、单位上三角矩阵的乘积。对于对称矩阵，可以分解成矩阵分解的形式，特别是对称正定矩阵，可以分解成形式（称为Cholesky分解），其中为下三角矩阵由Doolittle分解的唯一性可知这些形式的分解也是唯一的。这些不同的分解都可以通过矩阵乘积的方法取得。下面例2.3以对称正定矩阵为例，说明通过矩阵乘积的方法取得矩阵分解的不同形式。当然，当矩阵的Doolittle分解存在唯一时，这些不同的分解分别时唯一的，因此可以通过待定系数的方法取得，也即通过直接分解的方法取得这些分

21、解。例2-3： 由此可知：进一步可以考察矩阵左上方各阶顺序主子式与三角分解所得矩阵的左上方的各阶顺序主子式的关系：若记矩阵A的各阶顺序主子矩阵为，同样的记号也适用于矩阵，则有，从而矩阵的各阶顺序主子式的值等于的相应的顺序主子式的值：，（因为）。以例2-3矩阵分解为例，容易得到：由此，也很容易看到，一个对称矩阵通过消去过程得到的分解（其中为单位下三角，为上三角矩阵），当的对角元全部为正数时，该对称矩阵必为对称正定矩阵。由消去的过程可知各次主元是（矩阵的对角元），可知当矩阵的各阶顺序主子式均非零时，（原始的）消去法可以顺利完成，这是因为当各阶顺序主子式均非零时，均非零，即消去法可以顺利完成。因此得

22、：定理2.1 若阶方阵的各阶顺序主子行列式均非零，则存在唯一的分解，其中为单位下三角矩阵，为上三角矩阵。进一步，当矩阵非奇时，列主元消去法必可顺利完成：因为当矩阵非奇时，的第一列必不全为零，故通过选主元，可进行第一步消元，即算法可执行，得到，且其下方全为零，因为所以的代数余子式也非零，（），即矩阵非奇，因此下一步列主元消去可进行，由此，可完成全部消元过程。2.2.4解三对角矩阵的追赶法定理2.2 若分解), 则(下)带宽,(上)带宽.证明: 对二和三阶矩阵显然.(确定),对矩阵的阶作归纳证明:设对阶矩阵结论成立.令,可验证(介绍分块运算):由归纳假设, 因此最后的矩阵乘法：前者是初等矩阵左乘实

23、际上是Gauss变换矩阵相乘，后者按分块运算容易得到。特殊情况,三对角矩阵：， 由前述的方法，很容易将它分解成两个特殊的下、上三角矩阵的乘积： ，其中因为未讲矩阵的直接分解，此处可由确定分解形式、待定系数方式取得：从而，在解方程组，其中时，可以将该方程组转化为求解两个简单方程组： 通常被称为“追赶法”： 。这只是Gauss消去法的一个具体应用，在计算机上的应用可以节约存储空间，减少运算量。若手工计算，实际只需应用Gauss消去法即可。例2.3 线性方程组解的可靠性向量与矩阵范数，误差的代数表征3.1向量与矩阵范数向量范数由二维或三维向量的长度概念推广而来工具。1、向量范数 向量，与之相应的一个

24、非负实数，记作，如果它满足条件：1） 非负性 ，2） 正齐性 ，3） 三角不等式 ；常用的范数1-范数 ，2-范数 ，-范数 ， 注：一般若不指某特定的范数，则范数符号不作下标，记为；2、矩阵范数 矩阵 ，与之相应的一个非负实数，记作，如果它满足条件：1） 非负性 ，2） 正齐性 ，3） 三角不等式 ，4） ，因为矩阵与向量相乘仍是向量，因此：特别要求一种矩阵范数与一种向量范数：5）矩阵与范数的相容性；与前述三种向量范数分别相容的常用的三种矩阵范数：1-范数 ，（最大列和）2-范数 ，-范数 ，（最大行和）；注：作为矩阵的特例向量，这些范数定义无矛盾。例：试证明：（a）若是范数，是一个非奇异矩

25、阵，则由定义了向量的一种范数；（b）对于矩阵，定义了与向量范数相容的矩阵范数。证明：a. 由于是范数，它必满足范数的三条件；由于，所以 非负性： 且 当且仅当 ，又由的非奇性，当且仅当时才有，因此：当且仅当； 正齐性： 三角不等式：因此，按此定义的范数是范数； b. 仿前，容易证明定义了一种矩阵范数。关于相容性：3、范数的等价性：4、矩阵谱半径： ；有(对任何一种范数) , 残向量与误差的代数表征（矩阵的条件数）：若记方程组(2.1)的准确解，近似解；则有近似解的误差：，近似解的残向量：。小小？例： 比较：方程组 与 ，从前者到后者仅是右端产生误差，而对应的解则由变为，即也产生误差；它们之间有

26、关系： ； 注意以上公式：由，以及 得到。此处，从本质上反映了方程组右端的相对误差导致解的相对误差的数值关系，反映了方程组原始数据发生扰动时，方程组解发生变化的大小。比较前已提到的“问题的性态”，数反映了方程组的性态。因此称此数为方程组或矩阵的“条件数”，记作 ；由于总有（实际上，除了外，几乎总有），而从推导可知，（此处“”表示“相当于”），因此一般总有； 上例中：对于更一般的矩阵和右端向量均发生扰动和的情况，有证明：存在，且 反证：若 奇异，即，则 有非 解 又矛盾，所以 非奇，存在，由 （至少对常用的三种范数）：由此得 .若: Hilbert Matrix nn3748596102.4线性

27、方程组的迭代解法迭代法：将方程组 （2.1）转化为等价的形式从而，将向量代入(2.4)的右边，可计算得到一个新的值：如此反复地进行，便得到向量序列；问题 I）如何将方程组（2.1）的形式转化为(2.3)？ II）给定初始向量，迭代产生的向量序列是否收敛？ III）向量序列的收敛性与初始向量的选取是否相关？ IV）如果收敛，如何估计误差？基本迭代法首先我们回答问题 I）:由方程组（2.1）：，使等式左端仅保留向量，其他一概放到右端，即：将方程组的第个方程（设，） 改成：，再将此式遍除便得：从而原方程组可改写成如下等价形式：Jacobi迭代：将代入上式右端，便可（按顺序逐行）进行计算得到 ，这便是

28、迭代：Gauss-Seidel迭代：按通常的串行计算方式，迭代（2.5）中先计算第一式得到,此数完全可以参与第二式的右端的计算,依次类推，便得到： 这就是迭代。2.4.2 迭代法的矩阵表示将方程组（2.1）的系数矩阵A分裂成三部分：，其中原始方程组形成，即左侧仅有若，分别对第方程除以，并记左侧上标为k+1，右侧上标为k，则有Jacobi迭代：若由“最新信息优先”，即将已生成的立即用于计算，便是Gauss-Seidel迭代：矩阵分裂：Jacobi迭代：Gauss-Seidel迭代：显然，不同的迭代方法是对矩阵A做不同的分裂得到的：Jacobi迭代：Gauss-Seidel迭代： 收敛性以下我们回

29、答问题II）和III）：记方程组的解为，若当时向量序列收敛于方程的解，即，有等价关系（注意到范数的等价性）：由于满足关系式 ，将它与一般迭代：，两者相减，可得 ，遍取向量范数，并由矩阵范数与向量范数的相容性，可得：从而，由此可见，当时，无论如何取，总有，即：；即定理2.3：迭代 取任意初值生成的向量序列当时收敛于方程的解的充分条件是.具体化考察对原始问题（2.1）进行迭代或迭代时，迭代收敛对矩阵的要求。由式可知，迭代的迭代矩阵为，迭代的迭代矩阵为 ；由于的第行是：（）；因此，显而易见当为严格对角占优时该行各元素的绝对值之和小于1，所以：；结论是：推论2.4： 若方程组（2.1）的系数矩阵严格对

30、角占优,则任取初值的迭代收敛.其他结论：定理2.5： 若方程组（2.1）的系数矩阵严格对角占优,则任取初值的迭代收敛.证：记 严格对角占优 ，由的有限性，记 由 及 相减即 设 则所以 由 ，得 证完定理2.6： 若方程组（2.1）的系数矩阵对称正定,则任取初值的迭代收敛；当正定时，迭代收敛，否则就不收敛。最后，有：定理2.8：迭代收敛的充分必要条件是矩阵的谱半径.例：关于充分性：Jocabi：令则因此，总有 ，但并不满足定理的充分条件。例：，解此方程的Jacobi迭代与Gauss-Seidel迭代有相同的敛散性；例：讨论下述方程的Jacobi迭代与Gauss-Seidel迭代的敛散性，解又，

31、 等价于解： ； 结论：Jacobi 迭代收敛，Gauss-Seidel迭代不收敛。结论：Jacobi 迭代不收敛，Gauss-Seidel迭代收敛。例： 正定，不正定， Seidel收敛，Jacobi不收敛。例：例：证明：时Seidel迭代收敛，仅当 时Jacobi 迭代收敛。证：i. 时正定：，综上，时正定，所以Gauss-Seidel迭代收敛；ii所以，（当且）仅当 时Jacobi 迭代收敛。 正定： 正定；综合正定条件，得结论。但这是在正定的前提下的结论，如果不正定，条件是否仍如此？例：，以表示解的Jacobi、Gauss-Seidel迭代收敛的充分必要条件；I：Jacobi 迭代：J

32、acobi 迭代收敛； II：Gauss-Seidel迭代：等价于解：Gauss-Seidel迭代收敛。 迭代终止的判据定理2.9: 迭代,若,则任意取初值生成的向量序列当时收敛于方程的解,且有误差估计: 或 ；证：由前定理2.6可知，迭代对任选初值总收敛，又由及 令，便得结论。注1：前者称为“先验估计”，后者称为“事后估计”；但由证明的过程可知，这些估计实际上是“误差界”的估计；迭代法通常要求计算获得的解满足给定的误差界： ：l 先验估计：对要求的误差界，由估计式 可计算得 ，从而，可将由此所得的作为迭代步数的最大控制数，取；l 事后估计：对要求的误差界，由估计式 可知，在迭代计算过程中，应

33、验证不等式是否满足，如果满足，迭代就可终止，取。在实际计算中，通常将上述两种终止判据同时使用。注2：对于系数矩阵严格对角占优的迭代，可用定理2.8定义的作为收敛因子（即取），当然此时应取-范数；注3：若已肯定迭代收敛，但收敛因子q难以获得，可根据：对适当大的确定q的估计值。第三章 数据近似已知函数关系的数据点：，求一个适当的函数，它是已知函数族称为基函数的线性组合，即：寻求系数，使，与按某种标准相近似；通常按下述标准之一，使之取得极小：若取向量形式：,则上述三数分别为：；3.1 多项式插值一般，若有个数据点，上述近似标准要求即：，则称近似函数为满足插值条件(3.2)的插值函数，而点称为插值节点

34、。若满足插值条件(3.2)的函数的基函数取 ，则是多项式称满足插值条件(3.2)的多项式为插值多项式。代数学基本定理：n次多项式有且仅有n个根（包括重根）3.1.1 插值多项式根据插值条件(3.2)，可形成有个未知量的个方程：，这个方程组的系数矩阵称为Verdermonde矩阵，由其性质可知，当 时，该方程组的解存在且唯一。因此有定理3.1：对数据点，存在满足条件(3.3)的唯一的插值多项式；强调：插值多项式唯一 形式可以不同 Lagrange(形式)插值多项式 由两点直线方程谈起，；若有不超过次的多项式： ，可形成如下表：1000000001000000001000000001称为Lagra

35、nge基本插值多项式，显然它有个根：；由于次多项式有且仅有个根，因此注意到，可知：；若记 便有 不难证明 为所求的插值多项式，它被称为Lagrange（形式）插值多项式。缺点：增加插值节点后，Lagrange插值多项式发生大变化。3.1.3 Newton(形式)插值多项式定义：差商一阶差商： ；二阶差商： ；n阶差商： ；考虑逐次增加插值节点后，插值多项式的变化：记满足插值条件的（不超过）次插值多项式为，并记，即，注意：是次多项式；由故有 ；由于次多项式有且仅有个根，因此；确定系数：同样的方法可得：由此形成的插值多项式称作Newton插值多项式，记作：差商-性质1、对称性：对的任意排列，有证明

36、：前方法的Newton插值多项式的导出过程是：从仅满足一个插值条件的插值多项式开始，逐步增加插值条件，直到最后的，形成的，它的次项的系数就是；同样，若从仅满足一个插值条件的插值多项式开始，逐步增加插值条件，直到最后的，形成的的最高次（次）项的系数就必定是。由插值多项式的唯一性可知，这两个插值多项式是同一的，从而它们的次项的系数也是同一的，由此便得结论。 _例：求满足下述插值条件的插值多项式，解：差商表 -差商-性质2、设有直到n阶导数，则证明：设为插值节点，则对应的Newton插值多项式，令，显然，它有个零点：；由定理可知有n个零点，同理可知有n-1个零点，依此不难推出在中至少有一个零点，记为

37、，即：，而，由此得结论。例：多项式，且与点的选择无关；例：设 带导数条件的插值多项式在插值问题中有大量的问题是寻求满足导数条件的插值多项式，为此定义重节点差商设函数在处可导，则定义重节点差商(3.9)这是由导数的定义，有，上述定义也可以从差商性质2推出；由于令便得上述结论同理，由差商性质2还可推出，当有阶导数时 (3.10)注：由重节点差商定义可知，Newton插值多项式是Taylor展开式的推广:若在Newton插值多项式(3.9)中令,利用重节点差商便有这就是函数在点处的Taylor展开式.例：求满足下述插值条件的插值多项式，解：建立差商表：由此便可写出插值多项式 ；容易验证所得插值多项式

38、不仅满足函数值插值条件,而且满足导数值插值条件.3.1.5 插值公式的余项若将以为插值节点插值多项式，再增加一个插值节点，得插值多项式，由插值性质，有，据此，并由所取 的任意性，可将改作，便得插值公式的余项公式：利用差商性质2，及符号，又可以有余项的另一种表达形式例3.8：（第62页），以点为插值节点的插值多项式记为，求解：由余项公式：，令，便有差商表：因此： l Runge现象：问题：答1：若有界，有界， 则可。例：，显然 ，取区间，插值点等距分布：等分区间，步长 ，节点 逐步向外扩充，比较：使 最小：是区间中点，等分区间等分，则 而当或，则 所以， 从而 若插值点非等距分布，例如可证明：，

39、因此当 时，另一方面，即使 有界，当 无界时，仍可能无界。例：，显然，当经典的例子：函数，区间等分份，节点；（参考后图）注意：通常插值公式用到.附表：的函数插值多项式值0.840.860.880.900.920.950.962.67444.06913.94510.0471-10.3345-39.9524-50.86440.970.9730.9740.9750.9760.978-58.5447-59.5704-59.7279-59.7819-59.7254-58.2381解决方案：1) 分段多项式插值减小，低阶导数可估计、有界；同时保持小区间。2) 其他插值方法，例如用其他基函数插值等；3) 最

40、小二乘近似逼近点多，近似函数低阶（放弃插值）；4) 其他逼近方法。Runge 函数及插值的图形3.3 最小二乘近似 （线性）最小二乘问题的法方程已知数据点： ，求近似式：，使以下值取得极小：； (3.35)当时，该问题称为最小二乘问题，问题的解称为最小二乘解。通过对以为参数，以为极小化目标函数的计算，可以证明该最小二乘问题的解就是方程组的解，此方程组称为最小二乘问题的法方程组，其中 ；对此也可以从另一个角度作一简单的解释。作向量差 ；对于期望的最小值，当然最优是，但这也就是；由方程组理论可知，当方程组的方程个数大于未知量个数时（通常称为超定方程组），方程组一般无解，因此当时，方程组 一般无解。

41、但若将此方程组左乘，则方程组变成通常的阶方程组，从而可以求解，这个阶方程组就是法方程组。严格证明: 注意到 关于的极小与 对应极小是完全一样的.因此最小二乘问题的解应满足方程:为解释方便,将记为 ,注意到, 有交换求和顺序，便形成方程组：比较法方程 的具体表达式：可知：法方程的解 最小二乘问题的解。 正交化算法通常法方程组是病态的（），因此求解法方程组最好采用正交化方法，这是一个解线性方程组的稳定性非常好的方法。定义：若矩阵满足条件，即矩阵的逆矩阵是它的转置矩阵，则称此矩阵为正交矩阵。注意，从矩阵乘法规则容易看到，正交矩阵的列是互相正交的，而且它的列的2-范数长度为1。这也是它被称为正交矩阵的

42、原因。正交矩阵性质1：正交矩阵的乘积仍是正交矩阵。设是正交矩阵，则；正交矩阵性质2： （正交变换在二、三维空间中旋转变换坐标轴夹角总是直角不改变向量长度）；引理3.4：对任意非零向量，总存在正交矩阵，使，其中；证明：首先，对任意向量，若，容易验证矩阵为正交矩阵。其次，令 即 ；注意，向量是由其方向和长度唯一确定的，此向量 的方向是 ，而其长度则为：，因此取 。验证：，而注意到因此， ，这就说明了： ，即，对于给定的向量，取适当的向量，形成正交矩阵，使之满足引理。算法：记矩阵，针对矩阵的第一列，构造m阶正交矩阵,使：，从而有然后针对矩阵的第二列的后m-1个元素构造m-1阶正交矩阵，使之成为m-1

43、维空间的第一坐标方向向量：，再将扩充为m阶正交矩阵，以此左乘，就有照此办理，经过n步如此的正交化过程，就有注意到正交矩阵的乘积仍是正交矩阵的乘积，由此，法方程组可变化为：当的秩为n，即的n个列向量线性无关时，全不为0，对角矩阵满秩，它必有逆矩阵。因此，由，其中，两边左乘，得，解此方程组，便可得所求的；同时，根据正交矩阵的性质，可有 。由Lagrange插值多项式，待定函数方法构造余项：余项：设有直到阶连续导数,记 令 则 由Rolle定理，在包含全体 的区间中，有个零点，有个零点,有(至少)一个零点,记作,由于为不超过次多项式,因此例：求不超过2次插值多项式，使之满足条件：解：利用待定系数法，

44、令则由导数插值条件，应： ，立即解得：，因此所求.第四章 数值积分 ( Numerical Integration and Differentiation)记号： 积分， 求积公式，误差： ， ( Quadrature Formula )节点：， 求积系数：；4.1 内插求积, Newton-Cotes公式利用插值多项式： 积分;例如,插值多项式取Lagrange形式 , 便有及此类通过在节点 处满足插值条件的插值多项式导出的求积公式称为内插求积公式。其中：求积系数，误差显然，当函数为不超过次的多项式时，必有.由此可见,若函数为不超过某次的多项式，便有求积公式的误差，即.定义：若求积公式对任何

45、不超过m次的多项式成立：，而对次多项式等式不再成立，则称该求积公式的代数精度为m (事实上,可对一般的函数与其相应的近似公式定义代数精度).代数精度何时没有误差。本节讨论的公式的节点为等距分布：称为Newton-Cotes 公式 Newton-Cotes 公式1）n=1，梯形(Trapezoidal)公式 ，误差：利用积分中值定理（参见附录I，定理I.6：积分中值定理）2）n=2 Simpson公式，误差 注：该误差公式不能简单地对仅用此三点的插值多项式的余项公式积分获得。3）n=4 Cotes 公式 复化求积公式(Composite Numerical Integration)梯形公式、Si

46、mpson公式的 优点：简单，容易估计误差； 缺点：大区间 误差大改进：利用积分的可加性 ，1）复化梯形公式误差：若连续，由连续函数的介值定理(见例I.2)，及，结论：代数精度 1 ；计算精度一般，若有误差，误差与步长的关联程度；（ 若，则 ）2）复化Simpson公式误差：若连续，与复化梯形公式误差的推导一样，有结论：代数精度 3 计算精度 ； 步长的选取变步长积分法通常，对数值积分有一定的误差要求：给定误差界，使；1) 复化梯形积分法，由误差公式设在变化不大，即将以上两式相除，即得由此，得梯形积分的误差估计： ；从而，由可计算的值 估计；由于这是由约等式导出的，因此通常取 估计数值积分是否

47、达到误差要求。2) Simpson积分公式，由误差公式设在变化不大，即将以上两式相除，即得 ；从而，由 ；3) 变步长 比较步长之比为2:1的复化梯形积分公式： ， ，可得（理解）:由此可见，容易通过计算步长对分以后两个不同的复化梯形积分公式计算值的差，估计复化梯形积分公式的误差，由此确定是否达到所要求的精度。 Romberg积分由梯形积分的误差估计： ，设想将误差部分“归还”给梯形积分，应可以得到积分的一个更精确的近似值。通过计算，有：= = ；注意，梯形积分公式 的精度是，而通过适当的组合，可得到Simpson公式，它的精度是，精度得到很大的提高；同时代数精度也得到提高。这种方法并没有增加

48、函数计算量，仅增加了极少的代数运算，就能大大地提高计算精度，这类方法通常被称为“外推法”。用类似的方法，将Simpson公式适当地组合，可得到精度为，代数精度为5的 Cotes公式；再将Cotes公式适当地组合，可得到精度为，代数精度为7的 Romberg公式： ；在计算时，可列出如下的表进行计算：实际计算时，通常当时，计算终止，其中是事先要求满足的误差界；如果用相对误差作为终止判断，则取。一般, Romberg 方法用到 为止，尽管理论上还可以继续进行： 以取得精度更高的公式。但我们可以从两方面衡量该公式的实际意义。1、将此公式改为：，由于，因此可将此式认作是对的修正，而修正值是；当较大时，

49、与已是积分的比较准确的近似值，因此这部分的修正对的进一步改进影响一般比较小；2、注意到公式的误差：，公式的误差：，公式的误差：，公式的误差：，依次类推可知获得的值的误差应是，而的性态是很难确定的（如同插值多项式余项估计），因此，一般不再继续进行此类外推（Romberg）方法，仅用到 为止。例 积分：例： 积分： （效果并不明显）由此，也可从此体会到，Romberg 方法的效果（准确性），与被积函数的高阶导数的性态有关。 待定系数法例：求确定以下公式的系数 使之具有尽可能高的代数精度，求其代数精度，并求其误差： 解：取插值多项式，先求差商表：由此，通过对该插值多项式积分，得；因此求积公式为：；误

50、差：该插值多项式的余项：；由于在积分区间中，函数非正，因此由中值定理显然，这是积分公式。而且从公式的导出过程可以看到，公式实际上隐含着函数在积分区间中点处的导数值信息，因此作为等距节点数值积分公式，它用到了函数的4个信息，它具有4阶精度（），代数精度达到3.注：对于一般区间的积分，可用变量代换 将积分变换到区间：，再利用所获得的结果，便可得到前获得的积分公式，及其误差公式。但这样求积分公式比较麻烦，以下介绍另一种方法待定系数法。通常情况下，公式的代数精度越高，该公式越精确。因此，可以利用代数精度刻画的积分与求积公式的等式关系，确定未知的求积公式中的某些待定的参数，使公式达到尽可能高的代数精度。

51、前面定义代数精度：若求积公式对任何不超过m次的多项式成立：，即，而对m+1次多项式则不再成立，则称该求积公式的代数精度为m. 由前已得的梯形、Simpson等公式的误差表达式（一般形式为）容易判定其代数精度分别为1和3. 然而如果仅有公式,一般不易判定其代数精度, 因为无法用无穷多的“任何不超过m次的多项式”去验证是否总有。为此我们可以给出代数精度的另一种定义：求积公式的代数精度为，当且仅当对任何不超过m次的单项式成立：，而对则（参见本章习题1）。显然，这种定义只须次验证即可完成，从而确定求积公式的代数精度。另一方面，上述等式关系形成了个方程，从而可以确定求积公式的个系数。这样导出公式的方法称为待定系数法。此外，若一个求积公式的代数精度为，由代数精度的新的定义方式可知，公式误差的一般形式必为，因此确定求积公式的误差实际上是确定系数，若令，注意到的阶导数为常数：，因此由： 容易得到系数（注意：和在验证该公式的代数精度时都已经计算过），从而便可得到误差表达式。例（重复前例）：确定以下公式的系数 使之具有尽可能高的代数精度，求其代数精度，及其误差： ；解：取分别代入公式，由代数精度的等价条件，可以形成4个方程，据此确定公式的4个系数：解此方程