计算方法课件_易大义主编

数值计算方法,数理学院数学系,2,诚信声明,本课件中大量采用网络及其他渠道搜集的相关文字信息和图片信息，这些信息因数量较大，无法一一列明出处，本人在此郑重声明：这些相关资料的版权归原作者所有，本文引用仅仅用于教学目的。如有不妥，请与本人联系，联系方式：xznuckj。在此对相关资料的作者所付出的辛勤劳动表示衷心的感谢，并对作者表示诚挚敬意！本人郑重承诺：尊重知识，尊重劳动，尊重版权，学术诚信。,3,联系方式,课程网站/start/StudyNA/(校内访问)Email地址xznuckj办公室泉山9701室,4,本课程成绩的组成,平时成绩（占15）：包括出勤、课堂提问、讨论情况等。实验成绩（占25）：包括出勤、实验报告（预习报告）等。期末成绩（占60）。总成绩：100。,5,课程基本信息,总学时56其中上课40学时；上机16学时。主要讲授如何解决各种工程技术问题，用数学语言描述问题，即建立数学模型，将之转化为一个数学问题，寻求合适的近似计算方法，编程计算，充分发挥计算机的记忆和快速运算功能，寻求最佳方案。本课程先修课程为：高等数学、线性代数、程序设计语言,6,教材,计算方法,易大义等,浙江大学出版社,2002年第2版,7,主要参考书,1.数值分析引论,易大义陈道琦,浙江大学出版社,1998年第1版,8,主要参考书,2.数值分析基础教程,李庆杨,高等教育出版社,2001年第1版,9,主要参考书,3.数值方法和MATLAB实现与应用,（美）GeraldRecktenwald著伍卫国万群张辉等译,机械工业出版社,2004年第1版其他各类有关“数值分析”和“计算方法”的书,10,计算方法课程体系,第一章数值计算中的误差第二章插值法第三章曲线拟合的最小二乘法第四章数值积分第五章非线性方程的数值解法第六章方程组的数值解法第七章常微分方程数值解法,11,计算方法课程体系,本课程的内容,数值逼近,数值代数,常微分方程的数值方法,插值法数据拟合的最小二乘法数值积分和数值微分*,线性方程组的求解非线性方程组的求解矩阵特征值*,第一章数值计算中的误差,3学时,13,本章内容,1.1引言1.2误差的种类及其来源1.3绝对误差和相对误差1.4有效数字及其与误差的关系1.5误差的传播与估计1.6选用算法应遵循的原则小结作业与实验,14,本章要求,1.熟悉计算方法在解决实际问题中所处的地位,熟悉计算方法是以计算机为工具求近似解的数值方法；2.熟悉绝对误差（限），相对误差（限）及有效数字概念；3.熟悉公式；4.熟悉选用算法应遵循的原则。,15,1.1引言,解决科学技术和工程问题的步骤:什么是数值计算方法:将所预求解的数学模型简化成一系列算术运算和逻辑运算,以便在计算机上求解,并对算法的稳定性、收敛性和误差进行分析。,实际问题数学问题提供计算方法程序设计上机计算结果分析,16,1.1引言,简单地说，就是研究如何用计算机有效地解决一个数学问题。如何理解这两个含义？,这句话有两个含义,（1）有一个有效的数学方法,（2）一个能实现方法的有效程序（算法）,先看两个例子,17,1.1引言,算法影响计算的速度和效率（见课本P2秦九韶算法）例1古代中国人的贡献多项式的计值：设f(x)=a0 xn+a1xn-1+an-1x+an原始的算法需：n+n-1+1=n(n+1)/2次乘法。秦九昭算法：f(x)=(.(a0 x+a1)x+an-1)x+an仅需n次乘法。计算代价快速下降。,18,1.1引言,算法影响计算的精度例2设多项式为(x-2)9,我们来计算其在区间1.92,2.08上的值。令p(x)=(x-2)9q(x)=x918x8+144x7672x6+2016x5-4032x4+5376x34608x2+2304x-512则p(x)=q(x)，以下我们分别作画p(x)与q(x)的图。,19,上例说明，即使数学上的恒等公式，用计算机来算，结果也是不一样的。,p(x),q(x),20,1.2误差的种类及其来源,一.误差来源,例1，例2的结果的根源,实际问题,数学模型,建立算法,上机计算,结果,(初值误差)观测误差模型误差,(方法误差)截断误差,舍入误差,21,1.2误差的种类及其来源,二.误差分类1.模型误差（描述误差）/*ModelingError*/简化，抽象问题后建立的数学模型与实际问题之差。2.观测误差/*MeasurementError*/观测和实验得到的参量（物理量为电压、电流、温度等）,22,1.2误差的种类及其来源,3.截断误差（方法误差）/*TruncationError*/有限过程代替无限过程的误差（无穷级数求和，只能取前面有限项求和来近似代替）。这种计算方法本身出现的误差，所以也称为方法误差。如右端是截断误差。,23,1.2误差的种类及其来源,4.舍入误差/*RoundoffError*/计算机字长有限，一般实数不能精确存储，于是产生舍入误差。例如：在10位十进制数限制下：舍入误差很小，本课程将研究它在运算过程中是否能有效控制。,24,1.2误差的种类及其来源,据说,美军1910年的一次部队的命令传递是这样的:营长对值班军官:明晚大约8点钟左右,哈雷彗星将可能在这个地区看到,这种彗星每隔76年才能看见一次。命令所有士兵着野战服在操场上集合,我将向他们解释这一罕见的现象。如果下雨的话,就在礼堂集合,我为他们放一部有关彗星的影片。值班军官对连长:根据营长的命令,明晚8点哈雷彗星将在操场上空出现。如果下雨的话,就让士兵穿着野战服列队前往礼堂,这一罕见的现象将在那里出现。,25,1.2误差的种类及其来源,连长对排长:根据营长的命令,明晚8点,非凡的哈雷彗星将身穿野战服在礼堂中出现。如果操场上下雨,营长将下达另一个命令,这种命令每隔76年才会出现一次。排长对班长:明晚8点,营长将带着哈雷彗星在礼堂中出现,这是每隔76年才有的事。如果下雨的话,营长将命令彗星穿上野战服到操场上去。班长对士兵:在明晚8点下雨的时候,著名的76岁哈雷将军将在营长的陪同下身着野战服,开着他那“彗星”牌汽车,经过操场前往礼堂。,26,1.3绝对误差和相对误差,一绝对误差/*absoluteerror*/设准确值，近似值。称为的绝对误差(简称误差)为的绝对误差限。二相对误差/*relativeerror*/称为的相对误差。实用中，常用表示的相对误差。称为的相对误差限。,27,1.4有效数字及其与误差的关系,一有效数字/*significantdigits*/,一定要从规格化后的数来判断其位数有效位数与第一个非0项后的数字个数是不一致的。四舍五入所得到的数是一致的。,28,1.4有效数字及其与误差的关系,注:四舍五入规则修正为四舍五以上入,五时若前一位是偶数则5舍去,奇数则进一,以减少积累误差。,29,1.4有效数字及其与误差的关系,30,1.4有效数字及其与误差的关系,二有效数位与误差的关系,31,1.4有效数字及其与误差的关系,证：,32,1.4有效数字及其与误差的关系,例5：为使*的相对误差小于0.001%,至少应取几位有效数字？解：假设*取到n位有效数字，则其相对误差上限为要保证其相对误差小于0.001%，只要保证其上限满足已知a1=3，则从以上不等式可解得n6log6，即n6，应取*=3.14159。,33,1.5误差的传播与估计,例：蝴蝶效应纽约的一只蝴蝶翅膀一拍，风和日丽的北京就刮起台风来了？！以上是一个病态问题/*ill-posedproblem*/关于本身是病态的问题，还是留给数学家去头痛吧！,34,1.5误差的传播与估计,一.一元函数情形,35,1.5误差的传播与估计,二.多元函数情形,36,1.5误差的传播与估计,二.多元函数情形（续）,37,1.5误差的传播与估计,例6：测得某桌面的长a的近似值a*=120cm,宽b的近似值b*=60cm。若已知|e(a*)|0.2cm,|e(b*)|0.1cm。试求近似面积s*=a*b*的绝对误差限与相对误差限。,38,1.6选用算法应遵循的原则,一.尽量简化计算步骤，减少乘除运算的次数。,39,1.6选用算法应遵循的原则,二.防止大数“吃掉”小数当|a|b|时，尽量避免a+b。例如，假设计算机只能存放10位尾数的十进制数，则例：,40,1.6选用算法应遵循的原则,算法1：利用求根公式在计算机内，109存为0.11010，1存为0.1101。做加法时，两加数的指数先向大指数对齐，再将浮点部分相加。即1的指数部分须变为1010，则：1=0.00000000011010，取单精度时就成为：109+1=0.100000001010+0.000000001010=0.100000001010,大数吃小数,41,1.6选用算法应遵循的原则,算法2：例：按从小到大、以及从大到小的顺序分别计算1+2+3+40+109,注：求和时从小到大相加,可使和的误差减小。,42,1.6选用算法应遵循的原则,三.尽量避免相近数相减例：a1=0.12345，a2=0.12346，各有5位有效数字。而a2-a1=0.00001,只剩下1位有效数字。,43,1.6选用算法应遵循的原则,四.避免绝对值很小的数做分母，会造成浮点溢出/*overflow*/当|b|a|时，应尽量避免。,44,1.6选用算法应遵循的原则,五.选用数值稳定性好的算法，以控制舍入误差高速增长,45,小结,1.1引言一.算法影响计算的速度和效率二.算法影响计算的精度1.2误差的种类及其来源一.误差来源二.误差分类1.模型误差（描述误差）2.观测误差3.截断误差（方法误差）4.舍入误差1.3绝对误差和相对误差一.绝对误差二.相对误差,46,小结,1.4有效数字及其与误差的关系一.有效数字二.有效数位与误差的关系1.5误差的传播与估计一.一元函数情形二.多元函数情形1.6选用算法应遵循的原则一.尽量简化计算步骤，减少乘除运算的次数。二.防止大数“吃掉”小数三.尽量避免相近数相减四.避免绝对值很小的数做分母，会造成浮点溢出五.选用数值稳定性好的算法，以控制舍入误差高速增长,47,作业与实验,作业（上作业本）：习题一（P26）：3、6、10实验实验名称：实验0截断误差与舍入误差实验目的：熟悉上机环境，了解截断误差、舍入误差的影响实验日期：自行练习，不单独安排实验时间具体要求另见文档,第二章插值法Interpolation,计算机学院陈克建5学时,49,本章内容,2.1引言2.2Lagrange插值多项式2.3Newton插值多项式2.补Hermite插值2.4分段低次插值2.5三次样条插值2.6数值微分（*）小结作业与实验,50,本章要求,1.熟悉插值法的含义及其几何意义；2.熟悉Lagrange插值公式及其余项的使用。3.熟悉差分的定义,会造差分表；4.会造差商表,并熟悉Newton插值公式的使用；5.熟悉差商与导数的关系式；6.熟悉简单的带导数条件的插值；7.熟悉分段插值法的含义。,51,2.1引言,本节内容一.插值问题提出二.几何意义三.插值多项式的存在唯一性返回章节目录,52,2.1引言,一.问题提出：表示两个变量x,y内在关系一般由函数式y=f(x)表达。但在实际问题中，有两种情况：1.由实验观测而得的一组离散数据(函数表),显然这种函数关系式y=f(x)存在且连续,但未知。2.函数解析表达式已知,但计算复杂,不便使用。通常也造函数表。如:y=sin(x),y=lg(x)。有时要求不在表上的函数值,怎么办？,53,2.1引言,P29,54,2.1引言,二.几何意义：两条曲线有交点（公共点）,P30,55,2.1引言,从几何上看,插值问题即:已知平面上n+1个不同的点(xi,yi)(i=0,1,2,.n),要寻找一条过这些点的多项式曲线。（不超过n次）问题：（1）满足插值条件的插值多项式P(x)是否存在？应该是几次多项式？（n次）（2）如果满足插值条件的多项式P(x)存在，应如何构造？（3）用插值多项式P(x)近似代替f(x)，误差如何？,56,2.1引言,三.插值多项式的存在唯一性,P31,57,线性代数“方阵的行列式”部分：线性方程组的系数行列式D0时，方程组有且仅有一个解,2.1引言,58,2.1引言,范德蒙(Vandermonde)行列式证明用数学归纳法证明,59,2.1引言,克莱姆(Cramer)法则克拉默法则如果方程组的系数行列式不等于零，即,60,2.2Lagrange插值多项式,本节内容一.线性插值二.抛物线插值三.Lagrange插值四.Lagrange插值多项式的存在性、唯一性五.截断误差六.算法返回章节目录,61,2.2Lagrange插值多项式,拉格朗日中值定理：,62,2.2Lagrange插值多项式,十八世纪法国数学家Lagrange对以往的插值算法进行研究与整理，提出了易于掌握和计算的统一公式，称为Lagrange插值公式。它的特例是线性插值公式和抛物线插值公式。线性插值抛物线插值Lagrange插值插值多项式的余项误差估计,63,2.2Lagrange插值多项式,一.线性插值已知两个插值点及其函数值：求一次多项式使得,插值节点,对应的函数值,64,2.2Lagrange插值多项式,由于方程组的系数行列式所以，按Gramer法则，有唯一解于是或（公式1）,P33线性插值公式,65,2.2Lagrange插值多项式,容易验证，过点(x0,f0)与(x1,f1)直线方程就是上式(公式1)，如下图所示。,66,2.2Lagrange插值多项式,例：已知解：用线性插值公式(公式1)，计算得到,67,2.2Lagrange插值多项式,应用条件：右图表明，对于象y=lnx这样连续光滑的曲线，当两个插值节点很近并且所求的函数值也很近时，用线性插值方法是足以保证精度的。,68,2.2Lagrange插值多项式,二.抛物线插值已知三个插值节点及其函数值：求二次多项式使得,插值节点,对应的函数值,69,2.2Lagrange插值多项式,由于方程组的系数行列式(3阶Vandermonde行列式)所以，有唯一解。即满足这样条件的二次多项式是唯一确定的。容易看出满足上述条件，所以它就是所求的二次多项式。,（公式2）P34抛物线插值公式,70,2.2Lagrange插值多项式,容易验证,P2(x)是过点(x0,f0)、(x1,f1)与(x2,f2)三点的抛物线，如下图所示。,71,2.2Lagrange插值多项式,例：已知解：用(公式2)，计算得到,72,2.2Lagrange插值多项式,应用条件：右图表明，对于象y=f(x)为连续光滑的曲线，当三个插值节点很近但所求的函数值却相差较大时，用抛物线插值方法是可以保证精度的。,73,2.2Lagrange插值多项式,三.Lagrange插值已知n+1个插值节点及其函数值：求次数不超过n的多项式Pn(x)。使得,插值节点,对应的函数值,74,2.2Lagrange插值多项式,由于方程组的系数行列式(n+1阶Vandermonde行列式)所以,这个n+1阶线性方程组,有唯一解,即Pn(x)是唯一确定的。容易验证:满足这些条件,所以Pn(x)就是所求的次数不超过n的插值多项式(存在性)。显然,式1、式2都是式3的特例。,公式3,75,2.2Lagrange插值多项式,插值多项式（公式3）插值基函数满足插值条件的多项式显然为：,P33,76,2.2Lagrange插值多项式,记则从而，对某一节点xk，有,77,2.2Lagrange插值多项式,Lagrange插值公式的标准型公式这样其中这改写了Lagrange插值公式，在许多理论分析中是非常有用的。,78,2.2Lagrange插值多项式,注：一次多项式插值过两点直线;二次多项式插值过三点抛物线;不用待定系数法(1)计算量大;(2)不易讨论误差;,79,2.2Lagrange插值多项式,四.Lagrange插值多项式的存在性、唯一性？首先我们回答两个问题：（1）满足条件式的插值多项式存在吗?（2）是唯一吗?定理满足插值条件式的插值多项式Pn(x)存在且唯一。,了解,80,2.2Lagrange插值多项式,证法1（1）存在性(采用构造性证明)构造基函数则,且令则显然有即我们找到了一个插值多项式Pn(x),81,2.2Lagrange插值多项式,证法1（2）唯一性(用反证法证明)假设有两个不同的插值多项式pn(x)和qn(x)存在令则且即有n+1个零点，所以与假设矛盾。证毕。,82,2.2Lagrange插值多项式,证法2,83,2.2Lagrange插值多项式,证毕。,84,2.2Lagrange插值多项式,五.截断误差：,P37Lagrange插值多项式的余项,85,2.2Lagrange插值多项式,P37,86,2.2Lagrange插值多项式,补充：,87,2.2Lagrange插值多项式,反复应用Rolle定理：,88,2.2Lagrange插值多项式,图应用Roll定理进行区间分析,89,2.2Lagrange插值多项式,90,2.2Lagrange插值多项式,在上例中,如果只给出两个节点169和225,也可以作插值多项式,即1次Lagrange插值多项式,有两个插值基函数,这种插值方法称为Lagrange线性插值,也可以在n+1个节点中取相邻的两个节点作线性插值。,91,2.2Lagrange插值多项式,Lagrange插值基函数为Lagrange线性插值多项式为,92,2.2Lagrange插值多项式,线性和二次比较，谁好？一般地，是否插值多项式的次数越高，其计算结果就越精确？,93,2.2Lagrange插值多项式,94,2.2Lagrange插值多项式,95,2.2Lagrange插值多项式,例4：教材P34例1P38例2P39例3,教材自学，此处不详解,自学，不详解,96,2.2Lagrange插值多项式,97,2.2Lagrange插值多项式,98,2.2Lagrange插值多项式,例6：,自学，不详解,99,2.2Lagrange插值多项式,100,2.2Lagrange插值多项式,Lagrange插值多项式的优点是：直观;对称;容易编程上机等。缺点是：插值基函数计算复杂;每增加一个节点,插值多项式的所有系数都得重算;计算上浪费。后面提出的Newton插值就是克服了以上缺点。,101,2.2Lagrange插值多项式,六.算法拉格朗日插值法N-S图：,i=1,i=n?,f=1,i=i+1,N,f=f*(X-xj)/(xi-xj),建立xn,yn,p=0,输入X,ji,Y,j=1,j=n?,j=j+1,p=p+f*yi,输出X和p,结束,102,2.2Lagrange插值多项式,Lagrange插值多项式求插值的Matlab程序,%lagrangen.mfunctiony=lagrangen(x0,y0,x)n=length(x0);m=length(x);fori=1:mz=x(i);s=0;fork=1:nL=1;forj=1:nifj=kL=L*(z-x0(j)/(x0(k)-x0(j);endends=s+L*y0(k);endy(i)=s;endy;,103,2.2Lagrange插值多项式,比较不同的插值多项式次数对插值的影响,%Chazhibijiao.mx=-5:0.1:5;z=0*x;y=1./(1+x.2);plot(x,z,k,x,y,r)axis(-55-1.52);pause,holdonforn=2:2:10 x0=linspace(-5,5,n+1);y0=1./(1+x0.2);x=-5:0.1:5;y1=lagrangen(x0,y0,x);plot(x,y1),pauseendy2=1./(1+x0.2);y=interp1(x0,y2,x);plot(x,y,k),holdoffgtext(n=2),gtext(n=4),gtext(n=6)gtext(n=8),gtext(n=10)gtext(f(x)=1/(1+x2),104,2.3Newton插值多项式,本节内容一.差分二.差商三.差商与差分关系四.Newton基本插值公式五.Newton插值计算步骤返回章节目录,105,2.3Newton插值多项式,Lagrange插值虽然易算，但若要增加一个节点时，全部基函数li(x)都需重新算过。将Ln(x)改写成的形式，希望每加一个节点时，只附加一项上去即可。下面首先讨论在实际问题中常遇到的等距节点(差分)情况，此时，牛顿插值公式被进一步简化。然后学习不等距节点(差商)的牛顿插值多项式。,106,2.3Newton插值多项式,一.差分,107,2.3Newton插值多项式,108,2.3Newton插值多项式,109,2.3Newton插值多项式,二差商（亦称均差）/*divideddifference*/差商是数值方法中的一个重要概念，它可以反映出列表函数的性质，并能对Lagrange插值公式给出新的表达形式，这就是Newton插值/*NewtonsInterpolation*/。,110,2.3Newton插值多项式,对已知的列表函数称为f关于xi的零阶差商，记为：由零阶差商出发，可归纳地定义各阶差商。,插值节点,对应的函数值,111,2.3Newton插值多项式,1.定义,一阶差商的差商,112,2.3Newton插值多项式,P47,113,2.4Newton插值多项式,3.差商表（实用）,规定函数值为零阶差商,P47,114,2.3Newton插值多项式,三.在等距节点的前提下,差商与差分有如下关系,115,2.3Newton插值多项式,116,2.3Newton插值多项式,依此类推,P41,117,2.3Newton插值多项式,4.例：对函数的列表部分，列出差商表（见表）f(x)=lgx差商计算表,118,2.3Newton插值多项式,四.Newton基本插值公式,牛顿公式,119,2.3Newton插值多项式,牛顿公式,P48,120,2.3Newton插值多项式,牛顿前差公式/*Newtonsforward-differenceformula*/牛顿后差公式/*Newtonsbackward-differenceformula*/,P42,P45,121,2.3Newton插值多项式,注：一般当x靠近x0时用前插，靠近xn时用后插，故两种公式亦称为表初公式和表末公式。,122,2.3Newton插值多项式,五.Newton插值计算步骤1.计算差商2.计算插值根据以上步骤,自己写出Newton插值算法流程图,并能参照流程图编程上机。,123,2.3Newton插值多项式,计算Newton插值步骤：首先根据节点和节点值计算差商表;然后利用Newton插值多项式估算f(x)。例：,124,2.3Newton插值多项式,解：先造差商表,125,2.3Newton插值多项式,由Newton公式得四次插值多项式为：例：教材P42例4，P45例5，P48例6,教材自学，此处不详解,126,2.补Hermite插值,本节内容一.问题提出二.数学描述三.求解思想返回章节目录,127,前述插值问题:要求被插函数与插值多项式在节点取相同值Lagrange型插值条件然而,不少实际问题不但要求在节点上函数值相等,而且还要求它的导数值也相等(即要求在节点上具有一阶光滑度),甚至要求高阶导数也相等,满足这种要求的插值多项式称埃尔米特(Hermite)插值多项式。,2.补Hermite插值,128,2.补Hermite插值,现代的仿生学就是一个典型的例子。在设计交通工具的外形，就是参照海豚的标本上已知点及已知点的导数，做插值在计算机上模拟海豚的外形制成飞机、汽车等外形。下面只讨论函数值与导数值个数相等的情况。,129,2.补Hermite插值,130,2.补Hermite插值,有关Hermite（埃尔米特）插值的详细方法等，此处略。具体请参阅教材或相关资料。,131,2.补Hermite插值,132,2.4分段插值,分段插值(piecewisepolynomialapproximation)本节内容一.Runge现象二.分段线性插值三.分段Hermite插值返回章节目录,133,2.4分段插值,一.计算中的Runge现象由插值问题的提出,通常我们会觉得当节点越来越密时,插值函数越来越接近于原函数。但是结果并非如此,因为多项式是上下震荡的,震荡的幅度不尽相同,不同区段的震荡密度也不一样,由此导致,利用较高阶的插值多项式所计算的结果,与原来的函数值相差甚远。这说明高次插值未必可行。结果表明,并不是插值多项式的次数越高,插值效果越好,精度也不一定是随次数的提高而升高,这种现象在上个世纪初由Runge发现,故称为Runge现象。,P50,134,2.4分段插值,不同次数的Lagrange插值多项式的比较图,Runge现象,135,2.4分段插值,例：,n越大，端点附近抖动越大，称为Runge现象,136,2.4分段插值,二.分段线性插值/*piecewiselinearinterpolation*/,137,2.4分段插值,三.分段Hermite插值/*Hermitepiecewisepolynomials*/,138,2.4分段插值,分段线性插值（图）,139,2.4分段插值,分段三次Hermite插值（图）,140,2.5三次样条插值,本节内容一.分段插值法二.三次样条插值三.三次样条插值函数四.边界条件五.三次样条插值函数的表达式六.例返回章节目录,141,2.5三次样条插值,一.分段插值法：问题：结点增多,多项式次数增高,逼近精度越好?未必!多结点高次插值往往在局部误差更大Runge（龙格）现象。实用：采用分段低次插值。有分段线性,分段二次插值等,其几何意义缺点：分段插值函数只能保证连续性,不能保证光滑性。,折线代替曲线,142,2.5三次样条插值,二.三次样条插值/*CubicSpline*/分段插值可以得到整体连续函数，但在连接点处一般不光滑，而Hermite插值虽然在连接点处一阶光滑，但整体插值由于结点多次数高而有可能发生龙格现象。希望：既想分段插值，又想在结点处保持光滑，甚至二阶光滑三次样条。,P53,143,2.5三次样条插值,什么是样条：是指飞机或轮船等的制造过程中为描绘出光滑的外形曲线（放样）所用的工具；样条本质上是一段一段的三次多项式拼合而成的曲线；在拼接处，不仅函数是连续的，且一阶和二阶导数也是连续的；1946年，Schoenberg将样条引入数学，即所谓的样条函数。,144,2.5三次样条插值,三.三次样条插值函数,P54,145,2.5三次样条插值,四.边界条件S(x)在区间xi-1,xi上是三次多项式,S(x)=aix3+bix2+cix+di,有4个待定系数,要确定S(x)共要4n个待定系数。由S(xi)=yi，i=0,1,n,有2n个条件。再由S(xi-0)=S(xi+0),i=1,2,n-1,有n-1个条件及S(xi-0)=S(xi+0),i=1,2,n-1,有n-1个条件共有4n-2个条件。尚缺少两个条件。为得到唯一的三次样条函数,通常可在区间a,b的端点x0=a,xn=b上各加一个条件,称为边界条件。,P54,146,2.5三次样条插值,注：4n-2个条件S(xi)=yi，i=0,1,n,有n+1个条件。S(xi-0)=S(xi+0),i=1,2,n-1,有n-1个条件S(xi-0)=S(xi+0),i=1,2,n-1,有n-1个条件S(xi-0)=S(xi+0),i=1,2,n-1,有n-1个条件,147,2.5三次样条插值,常用的边界条件有1.S(x0)=y0,S(xn)=yn;(夹持条件)2.S(x0)=y0,S(xn)=yn;(自然边界条件)3.假设(x)是以b-a为周期的周期函数,这时要求S(x0+0)=S(xn-0)S(x0+0)=S(xn-0)S(x0+0)=S(xn-0)这样确定的S(x)为周期样条函数。(周期性条件),P55,148,2.5三次样条插值,有关三次样条插值多项式的存在性和唯一性证明，此处略。通过具体求解待定系数来证明，具体请参考相关资料。,149,2.5三次样条插值,五.三次样条插值函数的表达式若假设S(xi)=mi,i=0,1,n,利用分段Hermite插值多项式,当xxi-1,xi时,有其中hi=xi-xi-1。为确定S(x),只需确定mi,i=0,1,n。可利用S(xi-0)=S(xi+0)来求出mi。,P55,150,2.5三次样条插值,151,2.5三次样条插值,所以,152,2.5三次样条插值,于是有由连续性条件S(xi-0)=S(xi+0)可得,153,2.5三次样条插值,3(ixi-1,xi+ixi,xi+1)=gi,则有imi-1+2mi+imi+1=gi,i=1,2,n-1.(*式)再结合不同的边界条件,可得关于mi的方程。,154,2.5三次样条插值,若边界条件为:m0=y0,mn=yn,带入(*式)可得,155,2.5三次样条插值,若边界条件为:S(x0)=y0,S(xn)=yn,则有连同(*式)一起,可得,156,2.5三次样条插值,若边界条件为周期性边界条件，由S(x0+0)=S(xn-0),和S(x0+0)=S(xn-0),有m0=mn和nmn-1+2mn+nm1=gn，其中于是有,157,2.5三次样条插值,对应不同的边界条件,只要求出相应的线性方程组的解,便得到三次样条函数在各区间xi-1,xi上的表达式由于三个方程组的系数矩阵都是严格对角占优矩阵,所以都有唯一解,前两个方程组均可用追赶法求解,第三个方程组可用LU分解法或Gauss消元法求解.,158,2.5三次样条插值,例：,先看教材P60例7,159,2.5三次样条插值,故有,160,2.5三次样条插值,例：,161,2.5三次样条插值,故有,162,小结,2.1引言2.2Lagrange插值多项式2.3Newton插值多项式2.补Hermite插值2.4分段低次插值2.5三次样条插值2.6数值微分（*）,163,作业与实验,作业（上作业本）：习题二（P68）：1、3、5实验实验名称：实验一拉格朗日插值法实验目的：验证拉格朗日插值多项式，进一步加深对拉格朗日插值法的理解。实验日期：09计11：2011年5月6日（周五）7、8节09计61：2011年5月3日（周二）1、2节具体要求见另外文档,第三章曲线拟合的最小二乘法,计算机学院陈克建4学时,165,本章内容,3.1引言3.2什么是最小二乘法3.3最小二乘解的求法3.4加权最小二乘法小结作业与实验,166,本章要求,1.熟悉插值法和拟合法的区别；2.了解偏差的概念;3.掌握使用最小二乘法进行数据拟合。,167,3.1引言,本节内容一.问题提出二.科学计算中两类逼近问题三.多项式逼近四.函数逼近问题描述五.插值和拟合的概念与区别返回章节目录,168,3.1引言,一.问题提出某种合成纤维的强度与其拉伸倍数有直接关系，下表是实际测定的24个纤维样品的强度与相应拉伸倍数的记录。提示：将拉伸倍数作为x，强度作为y，在座标纸上标出各点，可以发现什么?,169,3.1引言,170,3.1引言,从图中可以看出,纤维强度与拉伸倍数大致成线形关系,并且24个点大致分布在一条直线附近,可用一条直线来表示两者之间的关系。解：设y*=a+bxi,我们希望y*=a+bxi与所有的数据点(样本点)(xi,yi)越接近越好。即令=yi-y*i最小。必须找到一种度量标准来衡量什么曲线最接近所有数据点。,171,3.1引言,二.科学计算中两类逼近问题：1、关于数学函数的逼近问题：计算机只能做算术运算，因此，在计算机上计算数学函数必须用其它简单的函数来逼近，且用它来代替原来精确的数学函数的计算。如：f(x)=sin(x)用代替等。函数逼近的特点：（1）要求高精度逼近；（2）要求快速计算（计算量要小）。,无穷级数与函数逼近,172,3.1引言,2、建立实验数据的数学模型：给定函数的实验数据，需要用较简单和合适的函数来逼近（或拟合实验数据）例：已知y=f(x)实验数据希望建立y=f(x)数学模型（近似表达式）数据逼近的特点：（1）要求适度的精度；（2）实验数据有小的误差；（3）有些问题会有特殊信息来选择数学模型。,173,3.1引言,三.多项式逼近（已学过）1、Taylor多项式逼近函数（在xx0点）（详见教材P88）例：教材89例12、插值多项式逼近函数（详见教材P88，另教材第2章）,P88?,174,3.1引言,四.函数逼近问题描述设f(x)为a,b上连续函数,寻求一个近似函数P(x)(多项式)使在a,b上均匀逼近f(x)。,175,3.1引言,五.插值和拟合的概念与区别插值法是使用插值多项式来逼近未知或复杂函数的,它要求插值函数与被插函数在插值节点上函数值相同,而在其他点上没有要求。在非插值节点上有时函数值会相差很大。最佳逼近问题要求在被插函数的定义区间上,所选近似函数都能与被插函数有较好的近似。最佳逼近是在函数空间M中选P(x)满足,176,3.1引言,但由于绝对值函数不宜进行分析运算,常将上式化为来讨论,于是最佳逼近问题变为最佳平方逼近问题,而离散的最佳平方逼进问题就是常说的曲线拟合它们都可用最小二乘法求解。插值法适用于数据精确或可靠度较高的情况曲线拟合法适用于数据本身就有误差的情况,177,3.2什么是最小二乘法,本节内容一.问题背景二.偏差的概念三.最小二乘原则四.最小二乘法返回章节目录,178,3.2什么是最小二乘法,一.问题背景科学实验，统计分析，获得大量数据,179,3.2什么是最小二乘法,180,3.2什么是最小二乘法,当数据量特别大时一般不用插值法。这是因为数据量很大时所求插值曲线中的未知参数就很多，而且数据量很大时，多项式插值会出现高次插值（效果不理想）或分段低次插值（精度不高）；另外，测量数据本身往往就有误差，所以，使插值曲线刻意经过这些点也不必要。而曲线拟合是，首先根据物理规律或描点画草图确定一条用来拟合的函数曲线形式，也可选择低次多项式形式（所含参数比较少），然后按最小二乘法求出该曲线，它未必经过所有已知点，但它能反映出数据的基本趋势，且误差最小，效果比较好。,181,3.2什么是最小二乘法,二.偏差（残差）的概念,在回归分析中称为残差,仍然是已知x1xm;y1ym,求一个简单易算的近似函数y(x)f(x)。,但是m很大；,i本身是测量值,不准确,即if(xi),这时没必要取(xi)=yi,而要使(xi)yi总体上尽可能小。,P71,182,3.2什么是最小二乘法,常见做法：,使最小/*minimaxproblem*/偏差最大绝对值最小,使最小偏差绝对值之和,使最小/*Least-Squaresmethod*/偏差平方和最小,太复杂,不可导，求解困难,183,3.2什么是最小二乘法,三.最小二乘原则1.最小二乘原则使偏差平方和最小（上页中方法3）的原则称为最小二乘原则；2.最小二乘法按最小二乘原则选择拟合曲线y=(x)的方法,P71,184,3.2什么是最小二乘法,四.最小二乘法,P72,185,3.2什么是最小二乘法,186,3.3最小二乘解的求法,本节内容一.最小二乘解的求法二.最小二乘拟合多项式的存在唯一性三.一般最小二乘拟合返回章节目录,187,3.3最小二乘解的求法,一.最小二乘解的求法,188,3.3最小二乘解的求法,P733.3.2,189,3.3最小二乘解的求法,190,3.3最小二乘解的求法,P74m次多项式拟合,191,3.3最小二乘解的求法,（1）直线拟合（一次多项式拟合）若，a0,a1满足法方程组即a0,a1是法方程组的解。,192,3.3最小二乘解的求法,例1已知一组试验数据试用直线拟合这组数据.(计算过程保留3位小数)。解设直线ya0+a1x，那么a0,a1满足的法方程组公式为,193,故法方程组为解得a0=1.229a1=1.483所求直线方程为y=1.229+1.483x,3.3最小二乘解的求法,计算列表如下：,194,3.3最小二乘解的求法,（2）二次多项式拟合若满足法方程组即a0,a1,a2是法方程组的解。,195,3.3最小二乘解的求法,例2已知一组试验数据试用二次多项式拟合这组数据.(计算过程保留2位小数)解设满足的法方程组,196,3.3最小二乘解的求法,计算列表如下：,197,3.3最小二乘解的求法,故法方程组为解得a0=5.05a1=-4.04a2=1.01所求二次多项式为y=5.054.04x+1.01x2,198,3.3最小二乘解的求法,二.最小二乘拟合多项式的存在唯一性定理：设点x0,x1,xm互异，则法方程组存在唯一的一组解证明：用反证法（略）,199,3.3最小二乘解的求法,三.一般最小二乘拟合,200,3.3最小二乘解的求法,两种非线性最小二乘问题的求解途径（1）化为线性最小二乘问题部分可化为线性拟合问题的常见函数类型见下页表（2）马奎特（Marquardt）方法是在计算机上求解非线性最小二乘拟合问题的最为常用和有效的方法之一。（略）,201,3.3最小二乘解的求法,202,3.4加权最小二乘法,本节内容一.加权最小二乘法二.例返回章节目录,203,3.4加权最小二乘法,一.加权最小二乘法重度：即权重或者密度，统称为权系数。它的大小反映了数据（xi,yi）地位的强弱。定义加权平方误差为：,204,3.4加权最小二乘法,205,3.4加权最小二乘法,二.例例3已知一组实验数据（xi,yi）及权Wi如表所示。若x与y之间有线性关系ya+bx，试用最小二乘法确定系数a和b。,P84,206,3.4加权最小二乘法,解因拟合曲线为一次多项式曲线（直线）1(x)=a+bx，故相应的法方程组形式如式。将表中各已知数据代入即得法方程组解之得,207,小结,3.1引言3.2什么是最小二乘法3.3最小二乘解的求法3.4加权最小二乘法,208,作业与实验,作业（上作业本）：习题三（P89）：1、3实验实验名称：实验二最小二乘法实验目的：验证多项式曲线拟合的最小二乘解，进一步加深对最小二乘法的理解。实验日期：09计11：5月13日（周五）7、8节、20日（周五）7、8节09计61：5月10日（周二）1、2节、17日（周二）1、2节具体要求见另外文档,第四章数值积分,计算机学院6学时,210,本章内容,4.1引言4.2牛顿柯特斯公式4.3复合牛顿柯特斯公式（*）4.4龙贝格算法4.5高斯型求积公式（略）小结作业与实验,211,本章要求,1.理解数值求积的基本思想，掌握求积余项和代数精度的概念2.掌握梯形公式，Simpson公式及其误差估计；了解Cotes公式;3.掌握复合梯形公式，复合Simpson公式及其误差估计，了解复合Newton-Cotas公式；4.掌握龙贝格算法。,212,4.1引言,本节内容一.数值求积的必要性二.构造数值求积公式的基本方法三.求积公式的余项四.求积公式的代数精度返回章节目录,213,4.1引言,一.数值求积的必要性,214,4.1引言,P91,215,4.1引言,216,4.1引言,P91,217,4.1引言,P91,截断误差,218,4.1引言,P92,219,4.1引言,220,4.2牛顿柯特斯公式,本节内容一.牛顿柯特斯公式二.牛顿柯特斯公式余项三.牛顿柯特斯公式的数值稳定性和收敛性返回章节目录,221,4.2牛顿柯特斯公式,一.牛顿柯特斯(Newton-Cotes)公式,222,4.2牛顿柯特斯公式,P95,223,4.2牛顿柯特斯公式,P96,224,4.2牛顿柯特斯公式,P95,225,4.2牛顿柯特斯公式,226,4.2牛顿柯特斯公式,二.牛顿柯特斯(Newton-Cotes)公式余项,P97,227,4.2牛顿柯特斯公式,228,4.2牛顿柯特斯公式,三.牛顿柯特斯公式的数值稳定性和收敛性求积公式(4.1.3)的数值稳定性是指f(xk)的误差对数值积分结果的影响。若影响很大，则称该数值求积公式不稳定。求积公式的结点未必越多越好。因为有一个所谓的稳定性不能设到保证。另外，对于插值型公式而言，结点增多会因Runge现象而使求积误差增大。,229,4.2牛顿柯特斯公式,230,4.2牛顿柯特斯公式,231,4.2牛顿柯特斯公式,232,4.3复合牛顿柯特斯公式,本节内容一.复化数值求积法二.复化梯形公式三.复化Simpson公式四.复化Cotes公式五.误差估计六.复合求积公式步长的自动选取返回章节目录,233,4.3复合牛顿柯特斯公式,一.复化数值求积法提高求积精度增加节点分段使用节点少的Newton-Cotes公式即所谓的复化求积公式整体使用节点多的N-C公式。原因：高次插值有时出现Runge现象，误差更大；节点增多，Ak有正有负，不能保证稳定性。,P97,234,4.3复合牛顿柯特斯公式,例3：,235,4.3复合牛顿柯特斯公式,复化（复合）求积公式所谓复化求积，就是先将积分区间分成几个小区间，并在每个小区间上用低阶牛顿柯特斯公式计算积分的近似值，然后对这些近似值求和，从而得到所求积分的近似值。由此得到的一些具有更大实用价值的数值求积公式，统称为复化求积公式。,P97,236,4.3复合牛顿柯特斯公式,二.复化梯形公式,237,4.3复合牛顿柯特斯公式,P98,238,4.3复合牛顿柯特斯公式,三.复化Simpson公式,P98,239,4.3复合牛顿柯特斯公式,P98,240,4.3复合牛顿柯特斯公式,四.复化Cotes公式,P98,241,4.3复合牛顿柯特斯公式,例4：,P101,242,4.3复合牛顿