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

计算方法（C）第1章绪论数值计算数值方法的分析计算机上数的运算算法分析第2章线性代数方程组Gauss消去法消去法主元消去法矩阵分解Gauss消去法的矩阵意义矩阵的LU分解及其应用其他类型矩阵的分解解三对角矩阵的追赶法线性方程组解的可靠性向量和矩阵范数残向量与误差的代数表征解线性方程组解的迭代法基本迭代法迭代法的矩阵表示收敛性第3章数据近似多项式插值插值多项式Lagrange插值多项式Newton插值多项式带导数条件的插值多项式插值公式的余项2 最小二乘近似最小二乘问题的法方程正交化算法第4章数值微积分内插求积，Newton-Cotes公式Newton-Cotes公式复化求积公式步长的选取Romberg方法待定系数法数值微分插值公式方法Taylor公式方法（待定系数法）外推法第5章非线性方程求解解一元方程的迭代法简单迭代法Newton法割线法区间方法收敛性问题简单迭代一一不动点收敛性的改善Newton法的收敛性收敛速度第1章绪论数值计算现代科学的发展，已导致科学与技术的研究从定性前进到定量，尤其是现代数字计算机的出现及迅速发展，为复杂数学问题的定量研究与解决，提供了强有力的基础。通常我们面对的理论与技术问题，绝大多数都可以从其物理模型中抽象出数学模型，因此，求解这些数学模型已成为我们面临的重要任务。一、 本课程的任务：寻求解决各种数学问题的数值方法一一如何将高等数学的问题回归到初等数学（算术）的方法求解一一了解计算的基础方法，基本结构（否则只须知道数值软件）一一并研究其性质。立足点：面向数学一一解决数学问题面向计算机一一利用计算机作为工具充分发挥计算机的功能，设计算法，解决数学问题例如：迭代法、并行算法二、问题的类型1、离散问题：例如，求解线性方程组Ax=b——从离散数据：矩阵A和向量b，求解离散数据X；2、连续问题的离散化处理：例如，数值积分、数值微分、微分方程数值解;3、离散问题的连续化处理：例如，数据近似，统计分析计算；1.2数值方法的分析在本章中我们不具体讨论算法，首先讨论算法分析的基础一一误差。一般来讲，误差主要有两类、三种（对科学计算）：1）公式误差一一“截断误差”，数学一计算，算法形成一一主观（人为）：数学问题-数值方法的转换，用离散公式近似连续的数学函数进行计算时，一般都会发生误差，通常称之为“截断误差”；一一以后讨论2）舍入误差及输出入误差一一计算机，算法执行一一客观（机器）由于计算机的存储器、运算器的字长有限，在运算和存储中必然会发生最末若干位数字的舍入，形成舍入误差；在人机数据交换过程中，十进制数和二进制数的转换也会导致误差发生，这就是输入误差。这两种误差主要是由于计算机的字长有限，采用浮点数系所致。首先介绍浮点数系1.2.1计算机上的运算一浮点运算面向计算机设计的算法，则先要讨论在计算机上数的表示。科学记数法浮点数：约定尾数中小数点之前的数全为零，小数点后第一个数不能为零。目前，一般计算机都采用浮点数系，一个存储单元分成首数和尾数：XX士XXXXX首数I 尾数(t位)其中首数存放数的指数(或“阶”)部分，尾数存放有效数字。对于进制，尾数字长为t位的浮点数系F(P,t,L,U)中的(浮点)数，可以用以下形式表示：%d° d, . 1<%<Pfl(X)=±(H-+ )XPPp2 pt产 0<dj<P,j=2,3，…，t此处，指数l(称为阶)限制在L<l<U范围内。以下记实数系中的实数为xeR，它在浮点数系F(P,t,L,U)中对应的浮点数记为fl(x)eF(p,t,L,U)——p进制，t尾数位数，L,U阶的范围。几乎所有近代计算机都采用“二进制"(即P=2):位、字节和字分别是指位数不同的二进制数。例如十进制转换二进制11=200000000122=210000001044=220000010088=230000100099=23+20000010011010=23+21000010102727=16+8+2+1=24+23+21+2。00011011^ V ,1字节位是一个二进制数(即0或1)；字节是8个二进制数字；上表的最后一列是字节。单精度浮点数(singleprecision)按32位存储，双精度浮点数(doubleprecision)按64位存储。精度用于指明每个浮点数保留多少位以及尾数和阶数各分配多少位。单精度浮点数的尾数为23位、阶数为8位；双精度浮点数的尾

数为53位（包含符号位）、阶数为11位（包含符号位）。双精度浮点数的等价二进制数如下所示：64人位

f ^ ^bbbb…bb fddd…ddddd' V ' ^H' '11位指数（含符号位）符号位52位尾数浮点数的特点：1、实数转换到浮点数一浮点化，〈缺点：〉总会产生误差（除极个别的情况：x=21,l=0,±1,±2,…）按四舍五入，绝对误差：x-趴x）i<2p-（举例），〈优点:〉浮点化产生的相对误差有界（与数字本身的数量级无关）|3（x）|=x-fl（x）<1Pjx2注：设实数xeR，则按p进制可表达为：d、d、 dd 1<d<pTOC\o"1-5"\h\z±(1,2, ,t,t+1,)xPl 1x—±(7+^2 pt+产+…)*P 0<dj<P,j=2,3,…,t,t+1…按四舍五入的原则，当它进入浮点数系F(p,t,L,U)时，若d<1p，则t+1 2d1d° dfi(x)=±(4+-4+•••+-4)xpipp2 pt] dd、 d++1 ]若d>p，则fl(x)=±(+—++ +——)xpt+12 pp2 pt对第一种情况：, ,d 1 11cle|x-f(x)=(pt11+ )X0l< (—)xpl=p-pt+1 pt2 2|x一|x一fl(x)|=p—d 1+1—pt+1<7^(1)xpl=101T

pt2 2有|x|>-1pl，将此两者相除，便得p就是说总有：|x—fl(x有|x|>-1pl，将此两者相除，便得p另一方面，由浮点数要求的1<d1<p，x—fl(x) 1——-<-p1-1

x 22、每一个浮点数系的数字有限：2（p—1）pt—1（U—L+1）+13、浮点数系中的运算非自封闭，（因为数字有限、尾数字长有限、指数数字有限等）例：在F（10,4,-5,5）中，x=.5420x10-2,y=.2001x103，运算x*y和x/y,的结果显然已不在此浮点数系内，而x+y或x-y也不在此浮点数系内，需fl（x。y）结果才在此浮点数系内。浮点运算应注意：1）避免产生大结果的运算，尤其是避免小数作为除数参加运算；2）避免“大”“小”数相加减;3）避免相近数相减，防止大量有效数字损失；4）尽可能简化运算步骤，减少运算次数。原因：记A=maxB（x）|=maxx~f（x）=—P1-，由|x-fl（x）|<A|x|，可得：x2|（x。y）-f（x。y）|<A|x。y| （“。”表示任意一种四则运算）此处A是由机器字长（实质上是尾数字长/的大小）确定的常数（它反映了实际运算的精度）。显然，若要求运算的舍入误差小，应使运算结果（如：x土y,x.y,x/）较小。尤其是小分母运算：2-x=-，小yn大误差。yyy其次，当浮点数系中两个数量级相差较大的数相加（或减），注意到由于浮点数系中数字字长的有限性，可能导致“大数吃小数”。例如，F（10,4,-5,5）中x=.8231x103,y=.2317x10-1，则x土y=.8231x103±.2317x10-1=.8231x103±.0000x103=x似乎y（|y|<<|xI）没有参加运算。第三，同样，由于浮点数系中数字字长的有限性，当两个相近数相减时：例如，在F（10,8,-5,5）中，x=.82317844x103,y=.82317832x103，两数相减：x-y=.12000000x10-3，计算结果仅剩2位有效数字，而原来参加运算的数字有8位有效数字，这将严重影响最终计算结果的精度。1.2.2算法分析作为一个可用的算法，必须考虑其效率和可靠性，定义：计算机完成一个乘法和一个加法，称为一个浮点运算（记为flop）；

注：由于计算机在运算时，加(减)法所耗时间远少于乘(除)法，所以通常只须计算乘法的次数，因此也有“一个算法需要多少个‘乘法”这一提法。1、计算效率一可计算性(计算复杂性一空间、时间)例：解线性方程组Ax=b的Grammar方法：x,A4A「其中|A|是方程组系数矩阵A对应的行列式，而|A/则是以右端向量b替代A的第i列所得矩阵对应的行列式。由线性代数知识可知，若A=(aij)，则IA|二Z(-1)J(i；,i；…口aa…。.，1i2i '^li1 2 n其中J(i；,%…,in)是由｛1,2，…,n｝变换到｛i；,i2,…in｝所需的置换次数。可见每计算一个行列式，需要(n-1)-n!个浮点运算；因此，按Grammar方法解方程组约需N=(n2-1)-n!个浮点运算。当n=20时N氏9.70728x1020，用一个运算速度为108/秒的计算机进行求解，约需3.078x105年(日前报道我国计算机已达到3840x108/秒，这仍需近10年)。而n=20的方程组应该说是一个小型的方程组。因此Grammar方法是一个不能接受的算法，它缺乏可计算性。第二章将介绍的Gauss消去法和迭代法就有较高的效率，具有很好的可计算性。2、计算可靠性作为算法，除了考虑其效率外，必须重视可靠性，它包括两方面：问题的性态和方法的稳定性问题的性态所计算的问题当原始数据发生小扰动时，问题的解一般也发生扰动。问题的，性态——问题的解对原始数据发生变化的敏感性 。原始数据小扰动n原始数据小扰动n问题解小扰动一一问题是良态的大变化 问题是病态的例：线性方程组：TOC\o"1-5"\h\z\o"CurrentDocument"1 1 11x+—x+-x- 122336 m1 1 1 13 日V1\—x+—x+—x—— 的解是：X=1213243121 1 1 47 ⑴—x+—x+—x- 〔31425360若将方程组的系数改写成具有2位有效数字的小数:1.00x1.00x+0.50x+0.33x=1.8<0.50x1+0.33x2+0.25x=1.10.33x10.25x2+0.20x3=0.78i1 2 3(-6.22)的解则变成：X=38.25；

1-33.65J这是一个典型的病态方程组。一般：由原始数据Xn计算结果f(X),扰动后的数据~n计算结果f(~)，若对问题f存在常数m，满足关系式:f(x)-f(1)

f(X)f(x)-f(X)f(x)m=sup ~—x一~X则称(相对误差之比的上界)m为该问题的条件数，记作m=cond(f)；由微积分中值定理知识容易计算出，当xx~时，cond(f)exf^(X)。f(x)稍后我们在第二章将对线性方程组的性态作进一步的分析。算法的数值稳定性：例：计算I=e-1/xnexdx,n=0,1,,8；n0解：由微积分知识，可得计算公式，①I=1—nI，②I=i(1—I)，n n—1 n—1nn我们将准确值与计算结果列表如下：方法II1IIII £ III准确值0.63211.36792.26423.20734.17095.14556.1268 7 .11248.1010①.6321.3679.2642.2074.1704.1480.1120.2160-.7280②.6320.3680.2643.2073.1709.1455.1268.1124.1010由上表可见，方法①中，原始步的误差，随着计算步数的增加被严重地放大，特别是18竟变成负数(注意：被积函数是非负函数)，而方法②则相反；这是因为方法①中，若前步有误差6：~k1=1kl+3，则〜〜——I=1—kI=1—kI-k6=I—k6，说明Ik］的误差6，经一步计算后被扩大了k倍，随着k的增大，误差将被大大地扩大；而通过同样的分析可知：方法②中Ik的误差则被缩小k倍。算法的数值稳定性：算法对初始误差导致的最终误差的可控性 ，如果最终误差被有效地控制，则称算法是稳定的，否则就是不稳定的。

第二章线性方程组求解线性方程组：Ax=b%1Ax=b%1%1+%2%2a%+a%211 222+•••+a%1nn+•…+a%2nn=P1=P2(2.1)aan1%1+an2%2+…+an/n=pn其中--■12X--■12Xx:・x-X1222:n2aa♦a11•111CN«rga•a_=AGauss消去法消去法设计方法原则：面向计算机(事先未知元素，编制程序)，例:33n——2233n——232——%2基本思想：降维一N维问题转化为N-1维问题一逐次降维，依次进行消去过程一一对方程组(2.1)由上而下逐步消去对角元akk(k=1,2,…,n—1)以下的ak(i=k+1,…,n)，使之转化为如下等价形式，达到目标：a%+a%+…+a%=pa%++a%=p (22)…………^ J=Pn从而，可进入回代过程，再由下而上，逐步解得%k(k=n,n-1,…,2,1)：这儿的akk(k=1,2，…,n—1)——主元对问题(2.1)设a11H0,目标：将第2〜n方程的%1的系数a》，…,a雨转化为0；方法：“第k个方程”-%x“第1个方程”(k=2,3,…,n)，得到a11

'aj]+a12x2+…+a1x=01a(1)x+…+a(1)x=0(1)< 222 2nn2……………a(1)xh na(i)x=0(i)n22nnnn现在只须关心由第2〜n方程形成的n-1维方程组，以后可仿此进行。消去：第k步(k=1,2，…,n—1)：设a妹手0，以第k行为基准，消去以下各、- 一 a行中k列a以下的a(i=k+1,…,n),令l=一,施行：第k行一l义第i行kk ik ika kikkn新的第k行，元素变化：(a=0),a=a—la,0=0—l0,经过n—1kk ijijikkjiiikk步消去(注意：ann以下无元可消)，得到(2.2)式。〈注意：每计算1个鼠,a..,0i仅需1个浮点运算〉回代：第一步第二步x[0n-ax]n[0n-ax]n-1，nn/an-1,n-1[0k-(akk+1xk+1+・・・+aknQa'kkk=n-1,…,2,1;运算量:消元：n元方程组：第1步消元：对第i(i=2,3,…,n)行，共n-1行；每1行：计算li1(i=2,3,…,n),1个乘除法(或“浮点运算”)计算新的a(1)=a„-1,a.(j=2,3,,…,n);0⑴=0,-l*共(n-1)+1个TOC\o"1-5"\h\zj Ij111j I I11 1乘除法(或“浮点运算”)第1次元共(n-1)[1+(n-1)+1]=n2-1个乘除法，因此，消元的运算量£(k2-£(k2-1)=Z(k2-1)Wk2k=n k=1 k=1-n=—+————n3 2 6回代：第1步，求xn需要1个乘除法(或“浮点运算”)，即：第2步求xn-1用2个乘除法(或“浮点运算”)，一般，第k步求x「用k个浮点运算，因此，n-k+1回代的运算量才k=nn(n+1);k=1Gauss消去法的总运算量：T=竺+n2-n。3 3例2-1：解线性方程组

X1+2x2+x3=0<2X1+2x2+3x3=3—11—3x2解：例：k-1-3l=221—l=-131(1-1；解x=-2 1 3-1k1J(-2S(-2-22-13-1-1；解x=k-25Jk1-21Jk25-117-1k1-3-1181627641681256-2-3-2-112.1.2解：例：k-1-3l=221—l=-131(1-1；解x=-2 1 3-1k1J(-2S(-2-22-13-1-1；解x=k-25Jk1-21Jk25-117-1k1-3-1181627641681256-2-3-2-112.1.2404360J57J1044190J-147S-10-42J（列）主元消去法1JkS(11Jk-1S(1-1-3122424-2-3-1-11Jk10J1824J-20-74J-2-1(-1S；解x=；解x=-31J-1k1J算法中，若第k步：Z）a“kk或ii)|aJ«h,k\i=k+1,…,n则按照原来的简单Gauss消去法算法可能无法执行，也可能出现大误差:(10-5准确解:例：在浮点数系F(10,5,(10-5准确解:0.250001875x*=0.499998749解：按Gauss消去法解:(.10000*10-4.20000*1011.20000*101.30000*101.10000*101、 /21-020000*106~>.20000*101)(.10000*10-4 .20000*101 .10000*101( -.40000*106-.20000*106)X=(.00000•100.50000*100)T误差大，原因：若有误差ajnakj+5，Bk=Pk+5则(~=a-1a,p=p-1p，则演变为ijijikkjiiikka=a-1(a+5)=a+15,p=p-1(p+5)=p+15；ijijikkj ijik iiikk iik这说明前步各元素的误差均被放大1ik倍。克服方法，将按绝对值最大的元素交换到主元位置，使|1汝区1，使前步的误差不再被放大。(.10000*10-4.20000*101.10000*101)—^1行5行》1.20000*101 .30000*101.20000*101)(.20000*101 .30000*101.20000*101)———.50000*10-5~>(10000*10-4.20000*101.10000*101)(.20000*101.30000*101.20000*101)( .20000*101.10000*101)X=(.25000•100.50000*100)T消元过程中，第k步消元(k=1,2,…,n-1)以第k行为基准，消去其后的n-k个方程的a (系数矩阵第k列a以下的元素)，Gauss消去法要求主元ik kka卜产0。为避免出现akk=0作为主元，并使前步的误差不被放大，应使|1汝区1，为此应使：akkl=Max{aJak+1k|「・|aJ}，通常采用按列选主元的列主元消去法：若a"=Max{aJak+1J…卜J1便将第k行与第-行交换，使amk与akk交换位置:使新的第k行执行在原始Gauss消去法中的角色，保证将作为除数的主元akk|4|a,k|i=k+1,…，n，从而|1,k|41，重复前述的Gauss消去过程。列主元消去法的效果： ’

.算法稳定，即计算误差能被有效控制；.当矩阵A的行列式|A|丰0时，算法总可以执行完成；注：若矩阵A是对称正定或严格对角占优，则不选主元，消去法也是稳定的;矩阵严格对角占优的定义：对角元的绝对值大于该行其他元的绝对值之和，•>臼％尸1j丰i问题：为什么系数矩阵A严格对角占优时，不选主元的消去法也是稳定的?为保证消去法是稳定的，计算应如何进行?a—LjaTOC\o"1-5"\h\z\o"CurrentDocument"a aa—Lja注意：第一步消元6c=(~——出a=a-a—Lj，由于占优:jijaij ij iia11 11 11它的效果与主元消去法的I"<1一样，误差不会被放大。但因此算法也应适当改变，应先对第一行计算此值，然后从第2列起逐列从上到下计算。且第一步消元后生成的右下方(n-1)*(n-1)阶矩阵仍是严格对角占优，以第二列为例:+a211a+a211a—ja116cj=aj-a —j n11X—c12jj=3<X—2jj=3=X—2jj=312a11=XX—c12jj=3<X—2jj=3=X—2jj=312a11=X—12jj=3+X—11jj=3-a21又122aa-a —1211>a22-a21a—12a11a21/. |6c22vn>Xaj=3彳Xi1Jj=3a—12a11+xnj=3<〃j=3-a21+a21-a21a—12a11a—12a11即新的第二行的对角元绝对占优，其他也可同样证明。例2-2：列主元消去法求解例2-10、1一1121-2I-1I-13

%-/231

/2132、3 3%72%%3

%-/231

/2132、3 3%72%%2.2矩阵分解Gauss消去法的矩阵意义一一矩阵的三角分解若将Gauss消去法解方程过程中产生的lj（i>j）作为一个单位下三角矩阵——其对角元为1，对角线上的元素均为0——L的对应元素；将Gauss消去法消元过程最后得到的上三角矩阵一对角线以下元素均为0—记作U或U（仅考虑方程的系数矩阵部分）；如例2.1，再将L与U或U相乘:1 01/1 2LU=2-1LU=2-1、1/2112I-12-3LU=2 13J03J0显然LU相乘得到的正是原始所求解的方程的增广矩阵，而LU相乘得到的正是原始所求解的方程的系数矩阵。反过来，一个矩阵也可以通过Gauss求解的过程获得这样的“LU分解”，其中L为单位下三角矩阵，U为上三角矩阵。这样将一个矩阵分解成一个下三角矩阵和一个上三角矩阵的乘积形式，称为矩阵的三角（口00粒加）分解。1）消去法的矩阵意义： 以例2-1说明I'010I'010]02J210、-112J与以下矩阵乘法是一样的:210、-21--21-112)可见，一般的第一步消元是:L(A,b)=L(A,b)=1丫。11。11

aa11 11a11a11ano…o

=

、—1,

11

ppp,-1 1人aa…a' n1 八11 11 11若记第k步消元后形成的矩阵为Q（k）,b（k））,《特别：则第i步消元是将以下初等矩阵左乘矩阵（a（i-1）,b（i-1））形成31^31…~31L（a,b）=（A（1）,b（1））A（i）,b（。，1—l1i+1,k—lnk因此，Gauss消去过程从矩阵运算的角度来看，是LL^LL(A,b)=(U,y)其中U为上三角矩阵,y为向量：(%U产112R)产1n(y1]U=U22••巴2n:，y=[21^,nn；yn'矩阵分解A=LU注意到初等矩阵Lk具有性质:L-1

klni及L—L-1

klni及L—1L—1…L-12 2 n—11(11ln—1,2ln,2因此，我们有(A,b)=L—1L—1…L-1(U,y)=L(U,y)12 n—1根据矩阵乘积的性质，有A=LU这就是基本的矩阵三角分解一一Doolittle分解一一将矩阵分解为单位下三角矩阵L与上三角矩阵U的乘积。从矩阵乘法与行列式的关系可知，若矩阵A三角分解得到A=LU，则有：detA=detL*detU，由于下三角矩阵或上三角矩阵的行列式是全部对角元的乘积，因此可利用三角分解求矩阵对应的行列式的值。例如：上述例2.1方程组系数矩阵对应的行列式的值是-1，下例2.3对称正定矩阵对应的行列式的值为144。当系数矩阵为A的方程组可以顺利求解(不必选主元)时，解的过程显然是唯一的，因此它的Doolittle分解必然也是唯一的，从而可以通过待定系数的方法获得该矩阵的Doolittle分解(通常称为“直接分解”或“紧凑格式”)。其他三角分解除了上述的单位下三角矩阵与上三角矩阵的乘积形式以外，还有其他类型的分解，例如下三角矩阵和单位上三角矩阵的乘积，单位下三角矩阵与对角矩阵、

单位上三角矩阵的乘积LDMt。对于对称矩阵，可以分解成LDLt矩阵分解的形式，特别是对称正定矩阵，可以分解成GGT形式（称为Cholesky分解），其中G为下三角矩阵《由Doolittle分解的唯一性可知这些形式的分解也是唯一的》。这'4-20-2'4-20-22-30-313-7-2 012-31-31-1-71212-331491120-31312(411-2 1-330-3131(221212-331491120-31312(411-2 1-330-3131(2212-331110-31(2-110-3、2 3((2-102 1；2(313°724-20 4'(4-2 0 4'(4-20 4'(4-20l=-1。-2 2-3 121 21 -33l=-31-331-3-320-313-7l=0-313-7-42431l=3l=匕、4 1-7237l=141、 3-719,42、 210,43-2〈些不同的分解都可以通过矩阵乘积的方法取得。下面例2.3以对称正定矩阵为例，说明通过矩阵乘积的方法取得矩阵分解的不同形式。当然，当矩阵的Doolittle分解存在唯一时，这些不同的分解分别时唯一的，因此可以通过待定系数的方法取得，也即通过直接分解的方法取得这些分解。例2-3：由此可知:4(329)

进一步可以考察矩阵左上方各阶顺序主子式与三角分解所得矩阵的左上方的各阶顺序主子式的关系：若记矩阵A的各阶顺序主子矩阵为A(1),A(2)，…As)“，同样的记号也适用于矩阵L,U，则有A(k)=L^k)U(k)，从而矩阵A的各阶顺序主子式的值等于U的相应的顺序主子式的值：det(A(k))=det(U(k))(k=1,2,…,n)，(因为det(L(k))三1)。以例2-3矩阵分解为例，容易得到:4=1*4,(4 -2\o"CurrentDocument"-2 2、4=1*4,(4 -2\o"CurrentDocument"-2 2、0 -3(4[2(1(1(4-2-301-3-347由此，也很容易看到，一个对称矩阵通过Gauss消去过程得到的LU分解(其中L为单位下三角，U为上三角矩阵)，当U的对角元全部为正数时，该对称矩阵必为对称正定矩阵。由Gauss消去的过程可知各次主元是片],目22,…,目nn(矩阵U的对角元)，可知当矩阵A的各阶顺序主子式均非零时，(原始的)Guss消去法可以顺利完成，这是因为当各阶顺序主子式均非零时，%],目22,…,目nn均非零，即Gauss消去法可以顺利完成。因此得： ““定理2.1若n阶方阵的各阶顺序主子行列式均非零，则存在唯一的LU分解，其中L为单位下三角矩阵，U为上三角矩阵。进一步，当矩阵A非奇时，列主元Gauss消去法必可顺利完成：因为当矩阵A非奇时，A的第一列必不全为零，故通过选主元，可进行第一步消元，即算法L1P1A可执行，得至【」片产0，且其下方全为零，因为det(L1P1A)=detA丰0,所以％的代数余子式detA*也非零，(.「detA=det(L1P1A)=p11detA;产0),即矩阵A1非奇，因此下一步列主元Gauss消去可进行，由此，可完成全部消元过程。解三对角矩阵的追赶法上带宽q:j>i+q下带宽p:i>j+qXX・・・XX7定理2.2若A定理2.2若A=LU(Dooli分解)则L(下)带宽p,U(上)带宽q.证明：对二和三阶矩阵显然.（确定p,q），对矩阵的阶作归纳证明:设对n-1阶矩阵结论成立.令vp+vp+10…0),wT=(吗wq+1可验证（介绍分块运算）:'a(1

1—v1a1B vwawI'a(1

1—v1a1B vwawIn—1由归纳假设,B—1vwt

a=LU,因此(awBJf11

—vn-1人Y11v-VWawIn—1L1人wIn—1(11一vla最后的矩阵乘法：前者是初等矩阵左乘〈实际上是Gauss变换矩阵相乘〉，后者按分块运算容易得到。特殊情况,p=1,q=1——三对角矩阵:an-an-1bn-1n-1由前述的方法，很容易将它分解成两个特殊的下、上三角矩阵的乘积:TOC\o"1-5"\h\z(I C1T=LU=UT=LU=u1cn-1 n-1UJn其中〈因为未讲矩阵的直接分解，此处可由确定分解形式、待定系数方式取得〉:U1=b1,lkUk=bklkUk=bkkck-1,k=2,3,…，n从而，在解方程组Tx=d，其中d=（d1s2,…,盘产时，可以将该方程组转化为求解两个简单方程组：Tx=doLy=d,Ux=y,通常被称为“追赶法”：

y1=dy1=d1k=dk—lkyk-1k=2,3，…，nn/U< /n一一、/x=ykk kyk+1 k=n—1,…,2,1,k /Uk间，减少运算量。若手工计算，实际只需应用Gauss消去法即可。例27-5间，减少运算量。若手工计算，实际只需应用Gauss消去法即可。例27-52.3线性方程组解的可靠性——向量与矩阵范数，误差的代数表征向量与矩阵范数向量范数由二维或三维向量的长度概念推广而来一一工具。1、向量范数 向量x=(%,x2…乙/eRn，与之相应的一个非负实数，记作H，如果它满足条件：1)非负性1)非负性国＞0,||x||=0ox=0, VxeRn，2)3)正齐性桃||=卜|N，VaeR,xeRn，三角不等式||x+y||<||x||+||y||, Vx,yeR2)3)常用的范数1-范数2-范数口刈常用的范数1-范数2-范数口刈2=(k」2+以212+…+，1

xnl2)2，-范数注：一般若不指某特定的范数，则范数符号不作下标，记为|卜2、矩阵范数矩阵AeL(Rmxn)，与之相应的一个非负实数，记作I|a||2、矩阵范数如果它满足条件:1)非负性2)1)非负性2)正齐性||A||>0, ||A||=0oA=0, VAeL(Rmxn)，10a||=ai||a||,VaeR,VAeL(Rmxn)，3)4)三角不等式||A+b||<|All+IBII，VA,BeL(Rmxn)，IAbII^iaiiibii,VAeL(Rmxn),BeL(Rnp)，3)4)因为矩阵与向量相乘仍是向量，因此：特别要求一种矩阵范数与一种向量范数:5)|网引AIUHI,VAeL(Rm"),xeRn——矩阵与范数的相容性;

与前述三种向量范数分别相容的常用的三种矩阵范数:1-范数2-1-范数2-范数A1A|2=maxia芍，（最大列和）

ji=1__ 1-范数ML8ij（最大行和）;二（AtA的最大特征值）-范数ML8ij（最大行和）;j=1注：作为矩阵的特例一一向量，这些范数定义无矛盾。例：试证明：（a）若忡是范数，M是一个非奇异矩阵，则由何|=|吩||定义了向量xeRn的一种范数；（b）对于矩阵AeRn、n，|⑷=|M^M：||定义了与M向量范数11*11相容的矩阵范数。证明：a"由于卜||是范数，它必满足范数的三条件；由于||x||=||Mx『所以M⑴非负性：x=||Mr||>0,且||x||=||Mx||=0当且仅当Mx=0，又由M的非奇性，当且仅当x=0时才有Mx=0，因此：［x||=0当且仅当x=0；M⑵正齐性：怕||=im（ax）n=ia（Mx）n=aiiMxii=aiixiiM M⑶三角不等式：Ilx+y1M=|IM(x+y)||=1Mx+My\^11MxiI+1MyII=Ixl1M+IM1M因此，按此定义的范数x是范数;Mb.仿前，容易证明b.仿前，容易证明A1M=||Mam-1定义了一种矩阵范数。关于相容性:|MAM-iMx\<|MAM-“Mil=||A||||x||3、范数的等价性:Wl<M<^N，工WL<IWL<IWL，HI<IWI<MM，8 1 8pn1 2 1 8 2 84、矩阵谱半径：p(A)=max{\k|,九为A全体的特征值,i=1,2,-iii有(对任何一种范数) P(A)<||A||,...Ax=Xx,x丰0=|X|||x||<||Xx||=||A刈<||A“|x|n|九区Mgp(A)<||A||2.3.2残向量与误差的代数表征（矩阵的条件数）:若记方程组（2.1）的准确解x*，近似解~；则有近似解~的误差：e=x*-~=（ie2,…,e/，近似解~的残向量：r=b-A~。一小=向小？

例:2°00例:2°00Z°00〕x10.4991.001)'3.000、(1.500)"I(1)r=b-A~=(0.0015)比较：方程组Ax=b与Ax=b-r，从前者到后者仅是右端b产生误差-r，而对应的解则由x*变为~，即也产生误差e;它们之间有关系:群阳1A-1II1；1 ....1注意以上公式：由b=Ax=|b||=1Axi区lAlllxll=同<|A崎，以及r=b-Ax=Ax*-Ax=A(x*-~)=Aene=A-ir=回<M-i帆得到。此处，||a|||a-1从本质上反映了方程组右端的相对误差号导致解的相对误差J4的数值关系，反映了方程组原始数据发生扰动时，方程组解发生变化x*11的大小。比较前已提到的“问题的性态“，数M||at|反映了方程组的性态。因此称此数为方程组或矩阵的“条件数”，记作Cond(A)=|A|||a-证明：网<证明：网<1=G-BU存在，且%-B)-1||Cond(A)>1(实际上，除了A=I外，几乎总有Cond(A)>1)，而从推导可知，耳〜A||a-1|6

x*| "11耳〜A||a-1|6

x*| "11bli上例中：A-1=(此处“〜”表示“相当于”)，因此一般总有(1.001 -2.000、0.003(-0.4991.000)||A||||a-1||=3x3.001/0.003=30011< 11< 1-IBII反证：若G-B)奇异，即，则G-B)x=0有非0解对于更一般的矩阵A和右端向量b均发生扰动AA和Ab的情况，有Cond(A)< AA1-Cond(A)^AllnBx=x丰0又0=0=l(I一B)xW=llx-Bxll>IWI-llBxll>IIMHMx=(-网llxll>0矛盾，所以G-B)非奇，G-B)-1存在，由||I||=1(至少对常用的三种范数):=||I||= -B)-1(I-B)|= -BL-G-BLb|>M-B)-1||-|(I-B1||b||= -B)-11(1-||B||)由此得||G-b)-1||<金•

若：（A+AA）Q+Ax）=（b+Ab）=Ab-AA'x:.Ax=（A+AA）1Qb-AA.-x）

（A+AA）=A，+A-1AA）+A-1+A-1AA）1A-1（Ab-AA-x<——i-A-11kAA（Ab||+||AA||||x||）国＜M*阿f地w〕

M_巾啊［丽+wJ1bli=1同<1Alllxll .・・席V向闷<CondU） f幽W〕,llxll-1「4）1|AA||〔Ibll ||A||J1-Cond\A||A||HilbertMatrixHnCond（H）nCond（H）3 2 n 5.24x1027 2 n 4.75x10841.55x10481.53x101054.77x10594.93x101161.50x107101.60x10132.4线性方程组的迭代解法迭代法：将方程组Ax=b(2.1)转化为等价的形式x=Gx+d, (2.3)从而，将向量xk(k=0,1,…)代入(2.4)的右边，可计算得到一个新的值xk+1：xk+1=Gxk+d. (2.3‘)如此反复地进行，便得到向量序列Xkk,k=0,1,…)；问题I)如何将方程组(2.1)的形式转化为(2.3)？II)给定初始向%0，迭代(23)产生的向量序列^xk最否收敛？III)向量序列XkJ的收敛性与初始向量x0的选取是否相关？IV)如果收敛，如何估计误差？基本迭代法

首先我们回答问题I）:由方程组（2.1）：Ax=b，使等式左端仅保留向量x，其他一概放到右端，即：将方程组的第i个方程（设aii00,i=1,2,…,n）改成:ax=-ax-ax-iiii11.•-aax=-ax-ax-iiii11.•-aii-1xi-1-an+1x.+1再将此式遍除a..便得:iiax-ax-aiii11•.-aii-1xi-1-aii+1xi+1-ax+Pi)从而原方程组可改写成如下等价形式:A(a111(——(A(a111(——(-axa 21122-a12x2-a23x3-ax1nn-a2x+P1)+P2),(2.4)1(c C——\―ax-an22nn-1n-1a nn22nn-1n-1nnJacobi迭代：将xk=（xk,x2，…,xk）T代入上式右端，便可（按顺序逐行）进行计算得到xk+「（xk+1,xk+1,…,xk+W，这便是Jacobi迭代:1/xk+11/xk+1= (a111xk+1= (-axka21122-a12,•・-a】xk+Pp-axk+P2),(2.5)1xk1xk+1= (-axk-ana n11 n2nn-axknn-1n-1Gauss-Seidel迭代：按通常的串行计算方式，迭代Gauss-Seidel迭代：按通常的串行计算方式，迭代（2.5）中先计算第一式得到xk+1,此数完全可以参与第二式的右端的计算，依次类推，便得到:xk+1=—(

1a—axk—・・

xk+1=—(

1a—axk—・・

122・・一axk+P)1nnxk+1=<2工(一axk+1a22211—axk—・・・・一axk+P)2nn(2.6)1.xk+1=——(-axk+1-a

na n11nnxk+1—・・n22-axk+1nn-1n-1这就是Gauss-Seidel迭代。迭代法的矩阵表示将方程组（2.1）的系数矩阵A分裂成三部分：A=D-E-F，其中

/a )(0 0)(0a…a、1112 1naa00…aD=22••，-E=21...・.・・・・，-F=2n・・・・I aJ、 nn/[a]a2…0J[0 0J原始方程组Ax=b(ax+ax+…+ax=Pax+ax+…+ax=P《211 222 2nn 2 |ax+ax+…+ax=Pn11 n2U2 "n形成x=Gx+d，即左侧仅有xa”x1=-a[2x2 a[x+P]ax=—ax—ax —ax+PJ222 211 233 2nnr2… … … ax=-a1x1-a?x之 a】xI+P若a11a22…aj0，分别对第i方程除以a屋并记左侧上标为k+1,右侧上标为k,则有Jacobi迭代：xk+1 (axk a xk+P)1 a 122 1nn 1yxk+1一 (-axk-axk -axk+P)J2 a 211 233 2nnr2722… … 1xk+1— (——axk——axk——•———a xk+P)n a n11 n22 nn-1n-1 nI nn若由“最新信息优先”，即将已生成的xk+1l—1，2,…，i立即用于计算xk+1，便是l i+1Gauss-Seidel迭代：xk+1— ( -a,.kk a,xk+P,)1 a 122 1nn 111xk+1— (-axkk+1-axk-……-axk+P)J2 a 211 233 2nn 222… … … … … …1xk+1— (-axk+1-axk+1—・・・-a xk+1+P)n a n11 n22 nn-1n-1 nI nn0矩阵分裂：A—D-E-F("11(0 ]Ax-a22x2+a21x1. ,,lax ax+・・・+a x'nnn''n11 nn-1n-1/(a12x2+…讨1nxn]:+an-1nxn[ 0JJacobi迭代：Dx—Ex+Fx+bUx—D-1(E+F)x+D-1bUoxk+1—D-1(E+F)xk+D-1bGauss-Seidel迭代：(D-E)x—Fx+bUx—(D-E)-1Fx+(D-E)-1bUxk+1—(D-E)-1Fxk+(D-E)-1b

axk+1axk+1=—a111 12axk+1=—axk+1 —axk—•一-axk+B222 211 233 2nn2axk+1=—axk+1—a xk+1—・・・—a xk+1+Bnnn n11 n22 nn—1n—1 nxk+1=(D—E)—1Fxk+(D—E)―1bDxk+1—Exk+1=Fxk+b=Dxk+1=Exk+1+Fxk+b显然，不同的迭代方法是对矩阵A做不同的分裂得到的：Jacobi迭代：A=D—[E+F]nAx=boDx—[E+F]x=bnDx=[E+F]x+bnx=D-1[E+F]x+D-1bGauss-Seidel迭代：A=[D—E]—FnAx=bn[D—E]x—Fx=bn[D—E]x=Fx+bnx=[D—E]-1Fx+[D—E]-1b收敛性以下我们回答问题II)和III)：记方程组的解为x*，若当kT0时向量序列(力攵敛于方程的解x*，即limxk=x*，有等价关系(注意到范数的等价性)：ksxk-x*oxk—x*-0omaxxk—x*—0o卜—x*||-0o卜—x*1—0i 8由于x*满足关系式x*=Gx*+d，将它与一般迭代(24)：xk+1=Gxk+d，两者相减，可得xk+1—x*=G(xk—x*)，遍取向量范数，并由矩阵范数与向量范数的相容性，可得：卜+1-x*||<1GM-x*||从而，肾+1—x*<|G||||xk—x*<IGI|2|xk-1—x*|I<--|G||k+1卜0—x*||由此可见，当IGII<1时，无论x0如何取，总有卜k-x*|P0，即：xk了x*；即定理2.3：迭代xk+1=Gxk+d取任意初值x0生成的向量序列\kJ当kf8时收敛于方程的解x*的充分条件是gk1.具体化一一考察对原始问题(2.1)进行Jacobi迭代或Gauss-Seidel迭代时，迭代收敛对矩阵A的要求。由式(24)(25)可知，Jacobi迭代的迭代矩阵为D-1(E+F)，Gauss—Seidel迭代的迭代矩阵为(D—E)-1F；

由于D-1（由于D-1（E+F）的第i行是：（%aiiaii-1 0 0ii+1aiiaiia—in)；因aii此，显而易见当A为严格对角占优时该行各元素的绝对值之和小于1,所以:D-1（E+F）<1；结论是：推论2：4：若方程组（2.1）的系数矩阵A严格对角占优，则任取初值的Jacobi迭代收敛.其他结论:定理2.5：若方程组（2.1）的系数矩阵A严格对角占优，则任取初值的".号'j=i+1".号'j=i+1严格对角占优nl.+u,<1n—<1，由i的有限性,ii 1-i.i记r=maxui<1i1-li由 xk+1=D-1Exk+1+D-1Fxk+D-1b及x•=D-1Ex*+D-1Fx*+D-1b相减 ek+1=D-1Eek+1+D-1Fek即ek+1=-12(-a )ek+1+-1X(-a )eka ija ijiij=1 iij=i+1ek+i=maek+i=ek+iek+i=maek+i=ek+ik+1=ek+1<8 i0]艺it。业小i0i0j=11E|a||ek0j j001j=i0+1所以ek+1jXiai0i01j=1<-^0-8 1-li0ek+i叫00[4i0i01j=i0+1n00ekll00证完定理2.6：若方程组（2.1）的系数矩阵A对称正定，则i）任取初值的Gauss-Seidel迭代收敛；ii）当2D-A正定时，Jacobi迭代收敛，否则就不收敛。最后，有：定理2.8：迭代（23）收敛的充分必要条件是矩阵G的谱半径p（G）<1.例：关于充分性:♦+%2=例：关于充分性:X+10工=1l1 2Jocabi：xk+1=-1)xk||d-1(E+F)|=DT(E+F)|因此，总有例:la11、21。12]X=aJ22110-k1-10-k0但并不满足定理的充分条件。a11a22中0，解此方程的Jacobi迭代与Gauss-Seidel迭代有相同的敛散性;例：讨论下述方程Ax=b的Jacobi迭代与Gauss-Seidel迭代的敛散性，ii)22-1-12-10-1-2-20-22-10又，det(11-GS)=det11-(D-E11FLdetndet(入I-J)=入12解：]=12-211121211det11(D-E)-F=13-412+41=1(1-2)2np(GS)=2；结论：Jacobi迭代收敛，Gauss-Seidel迭代不收敛。ii) J=0-1 0/I/2-%、-10

JndetQI-J)=5 5. '=13+41np(J)=—detlx(D-E)-F]=2九2九-九-1

2九

-九12=813+812+21=21(21+1)212np(GS)=12结论：Jacobi迭代不收敛，Gauss-Seidel迭代收敛。例：A=3l4610正定，2D-A不正定，Seidel收敛，Jacobi不收敛。10例：A=4-3-2nlJll<1,13]15J)IJII<11例：A=证明：-1-2<a<1例：A=证：i.时Jacobi证：i.<a<1时A正定：A1=1>0,A3|=|A|=1-3a2+2a3=(1-a)2(1+2a)>0na>-%综上，<a<1时A正定，所以Gauss-Seidel迭代收敛;综上，-a)ii.-a-anii.-a-andet(Xl-J)=(k-a)2(1+2a)n p(J)=|2a|<1、-a-a所以，(1-a-a所以，(1-a-a1正定：det(2D-A)=-2a3-3a2+1=(a+1)2(1-2a)>0、一a-a（当且）仅当-L；<a<L；时Jacobi迭代收敛。.•.-1<a<1:正定；综合A正定条件，得结论。但这是在A正定的前提下的结论，如果A论，如果A不正定条件是否仍如此?(2a例：(2a例：A=b20a\KJVt-的充分必要条件；I：Jacobi迭代：以a,b表示解Ax=b的Jacobi、Gauss-Seidel迭代收敛J=DJ=D-1(E+F)=-%,加口1-J)=%k=k3—-abk4— 、；3/.p(J)— 、；3/.p(J)=——v|ab|n |ab|21Jacobi迭代收敛;II：Gauss-Seidel迭代：det(kI-GS)=detkI-(D-ELFLdetb0等价于解:2入 a2入 a0det11(D-E)-F]=b入2'b0a113=413—3ab12=412(1-jab)P(GS)=4abi|ab|<4oGauss-Seidel迭代收敛。迭代终止的判据定理2.9：迭代xk+1=Gxk+d，若冏=q<1,则任意取初值x0生成的向量序列Ik kT9时收敛于方程的解x*，且有误差估计：卜*-xklI-图卜1-x011 或k,-xklI-jzq卜+1-川|；证：由前定理2.6可知，迭代xk+1=Gxk+d对任选初值x0总收敛，又由VCm+1—xm—||Gxm—GxVCm+1—xm—||Gxm—Gxm-1—||G||||xm-xm-111—q||xm—xm-1-…qm-kxkk+1—xkxk1++1—xk<qkx1—x0|xj+1-xj11-(qm—k+qm—k-1+…q+1)||xk+1-xk11=1-qm-k+1xck+1-xk令m-8，便得结论。注1：前者称为“先验估计”，后者称为“事后估计”；但由证明的过程可知，这些估计实际上是“误差界”的估计；迭代法通常要求计算获得的解满足给定的误差界：卜*-引-£：•先验估计：对要求的误差界8，由估计式k ig((1-q)8x1-x0||)|x,-xk||<1q—|x1-x0||<8可计算得k> lgq-从而，可将由此所得的k作为迭代步数的最大控制数，取1—xk；•事后估计：对要求的误差界8，由估计式—”1<—”1<吉卜k+1-xk||<8 可知，在迭代计算过程中，应验证不等式卜k+1-xk||<8是否满足，如果满足，迭代就可终止，取~—xk。在实际计算中，通常将上述两种终止判据同时使用。注2：对于系数矩阵A严格对角占优的Gauss-Seidel迭代，可用定理2.8定义的R作为收敛因子(即取q=R)，当然此时应取8-范数；注3:若已肯定迭代收敛，但收敛因子q难以获得，可根据:-xm<qxm-xm-1对适当大的m确定q的估计值。第三章数据近似已知函数关系y=已知函数关系y=f(x)的数据点:它是已知函数族),IJ系数a『a?,…,a，使xi,y.)i=1m1求一个适当的函数p(x)，g｝(x)称为基函数一一的线性组合，即：寻求k(3.1)，P(x)=ag(x)+ag(x)+…+ag(x(3.1)，nn与f(x)按某种标准相近似；通常按下述标准之一，使之取得极小:E1=2|p(xi)-yi|i=1 E='Zi=1 E='Z(p(x)-y)22 iiii=14E=maxIp(x)-y8i=1,…,m若取向量形式:P=(P(x),P(x)，…，p(x))T,Y=(y,y,…，y)T,则上12 m述三数分别为:多项式插值(3.2)一般，若有n+1个数据点［xi,yi),i=0,1，…,n｝，上述近似标准要求p(x)=y,i=0,1,・・・，n(3.2)ii即：E1=E2=E8=0，则称近似函数p(x)为满足插值条件⑶2)的插值函数，而点｛(xi,yi),i=0,1,.,n｝称为插值节点。若满足插值条件⑶2)的函数p(x)的基函数取gk(x)=xk,k=0,1,2…｝，则p(x)是多项式(3.3)p(x)=a+axh naxn(3.3)称满足插值条件(3.2)的多项式p(x)为插值多项式。代数学基本定理：n次多项式有且仅有n个根(包括重根)插值多项式)的n+1个方程:根据插值条件(3.2)，可形成有n+1个未知量(a0,a「…,ap(x.)=a0+a1x.+…+anx.=y.i=)的n+1个方程:这个方程组的系数矩阵称为Verdermonde矩阵，由其性质可知,时，该方程组的解存在且唯一。因此有定理3.1对数据点1x.,y.),=0,1,....m，x,中xj,Vi中j｝，存在满足条件⑶3)

的唯一的插值多项式；强调：插值多项式唯 形式可以不同Lagrange(形式)插值多项式由两点直线方程谈起，%由两点直线方程谈起，%-Xc+ 0-y.%1-%0若有不超过n次的多项式/.(%)：l.(%.)=5,.=]0, i中J，可形成如下表:i 1JJ11, i=J%0%1…%i%n10(%)10(%)y01y000…000011(%)lj%)y1001y1…0000……………l.(%)l,(%)y.0000…1yi00………ln(%)i”(%)y„0000…001yn工l(%)yi iy0y1…yiyn§(%)称为Lagrange基本插值多项式，显然它有n个根：％0，…,%t-1，%1+1，…，%n；由于n次多项式有且仅有n个根，因此l(%)=a1(%-%)(%-%)…(%-%)(%-%)•••(%-%)注意到l.(%.)=1，可知：(%-%0＞一(％-%i-)(%-%i+J…(％-%”)(%1-%0)…(%1-%1-1)(%1-%1+1)…(%l-%n)'若记便有不难证明①(%若记便有不难证明①(%)=#Q-%.)=[ )% )...Q-%)i 0 1 nl.(%)1=03Q)二Q-%9&)L(%)=»l(%)y.i=0(3.4)(3.5)为所求的插值多项式，它被称为Lagrange(形式)插值多项式。缺点：增加插值节点后，Lagrange插值多项式发生大变化。

Newton(形式)插值多项式一阶差商:二阶差商:y[一阶差商:二阶差商:y[x0，x1,x2]二xi_x0y[x0,x2]y[x0,xj

x2-x1y[xn,x，…,x“Mx“]-y[xn,x，…,x“Mx“Jn阶差商：y[x,x，…,x]=0-1 n-2n01 n-2n-1-；01n xn-xn-1考虑逐次增加插值节点后，插值多项式的变化：记满足插值条件Lk(x)=y(i=0,1，…,k)的(不超过)k次插值多项式为Lk(x)，并记LJx)=L^1(x)+Q^(x)，即Qk(xQk(x)="x)-Lk-1(x)(Q0(x)=L0(x)=y0注意：Q(x)是k次多项式；由kLk(xi)=yi,i=0,1，…,k-1,k,Lk-1(xi)=yi, i=0,1，…,k-1,故有Qk(xi)=0,i=0,1，…,k-1；由于k次多项式Qk(x)有且仅有k个根，因此Qk(x)=ak(x-x0)(x-x1)^(x-xk-1)；确定系数ak：L0(x)=y0y0+ajx1-x0)=y1L1(x)=L0(xL0(x)=y0y0+ajx1-x0)=y1a=y——y0-=y[x,x],1x-x0110/.L(x)=y+y[x,x](x-x),类似L(x)=L(x)+Q(x)=y+y[x,x](x-x)+a(x-x)(x-x)=y,TOC\o"1-5"\h\z...a;.x0,x2]-y[x0,』=；x,x：x；「° 22021 2\o"CurrentDocument"2x-x 0122 1/.L(x)=y+y[x,x](x-x)+y[x,x,x](x-x)(x-x),同样的方法可得：Qk(x)=y[x0,%,…，xk](x-x0)(x-xj…(x-xk-1)由此形成的插值多项式称作Newton插值多项式，记作:N(x)=y+y[x,x](x-x)+y[x,x,x](x-x)(x-x)d—+y[x,x,…,x](x-x)(x-x)…(x-x) (3.6)

(3.7)差商-性质1、对称性：对瓜1，…,"的任意排列10,Z；…,J,有y[xo,",…，xj=y[x,o,x〃，…，x(3.7)证明：前方法的Newton插值多项式的导出过程是：从仅满足一个插值条件(xo,yo)的插值多项式Lo(x)=y0开始，逐步增加插值条件(\,yJ,(x2,y2)…，直到最后的(xn,yn)，形成的N(x)，它的n次项的系数就是y[x0,%,…,xn]；同样，若从仅满足一个插值条件(x,o,y,o)的插值多项式Lo(x)=y，o开始，逐步增加插值条件(x〃,y〃),(xi2,yi2)…，直到最后的(xn,y,n)，形成的N(x)的最高次(n次)项的系数就必定是y[x,1°,x〃,…,xin]o由插值多项式的唯一性可知，这两个插值多项式是同一的，从而它们的n次项的系数也是同一的，由此便得结论。x y(x) y[*,*]x y(x) y[*,*]y[*,*,*]y[*,*,*,*]xy\o"CurrentDocument"o oxy1 1xy2 2xy\o"CurrentDocument"3 o「ry一y

y[x,x]=-u o-o1x-x1oy-yy[x,x]=t p12x-xy2-y1y[x,x]=-3 232x—xy[x,x]—y[x,x]y[x,x,x]= 1——2 oo12 x-xy[x,x2]—y[x,x]y[x,x,x]= 2——3 1——3123 x—xy[xo,x1,x2,x3]y[xxxx]=y[],x2,x3]-y[x。,x、x2]

o,J2,3 x3-xo例：求满足下述插值条件的插值多项式，解：差商表nN(x)=-2+(x+1)+—(x+1)x-—(x+1)x(x-1)2 31 1-1-21 1-1-22 2o-1nN(x)=1+2(x-1)--(x-1)(x+1)--(x-1)(x+1)(x-2)

3 6 3差商-性质2、设y(x)有直到n阶导数，则乂x0乂x0,x1，…,x/二证明：设x0,x1，…，xn-1：—y(n)(己),m£[minx,maxx]n! iiii,xn为插值节点，则对应的Newton插值多项式(3.8)N(x)=y0+y[x0,x1Kx-x0)+ + y[x0,x]，…，xn](x-x0)(x-x1A…(x-xn_J，令R(x)=y(x)-N(x)，显然，它有n+1个零点：x0,x「…,xn1，xn；由Roll定理可知R'(x)有n个零点，同理可知R〃(x)有n-1个零点，依此不难推出R(n)(x)在[minx,maxx]中至少有一个零点，记为自，即:R(n)(己)二y(n)(己)-N(n)(己)二0，i i1而N(n乂9=n!y[x0,x^…,x]，由此得结论。例：多项式p(x)=Znaxknp[x,x，…,x]=a，且与点的选择无关;k=0k 0 1 nn例：设p(x)=x3+qx2+rx+s,ifp[0,1,t]=1,thenp[1,2,t]=?3.1.4带导数条件的插值多项式在插值问题中有大量的问题是寻求满足导数