插值与数值积分.ppt

大学数学实验大学数学实验Mathematical Experiments 实验实验3 3 插值与数值积分插值与数值积分1计算机会“算”吗？靠得住吗？例：把例：把4开开n次方，再平方次方，再平方n次，结果是次，结果是4？存在误差？存在误差？英国著名数值分析学家英国著名数值分析学家 HighamHigham(1998):(1998):Can Can you count on computers?you count on computers?精确计算：精确计算：解析结果解析结果(Analytical)近似近似计计算：算：数值数值结果结果(Numerical)n=55左右：结果变成左右：结果变成1计计算功效算功效=计算工具计算工具*计算方法计算方法(算法算法)浮点运算：舍入误差浮点运算：舍入误差2实验3的基本内容3.3.数值积分的数值积分的梯形公式、辛普森公式和梯形公式、辛普森公式和高斯公式。高斯公式。1.1.插值的基本原理；插值的基本原理；三种插值方法：拉格朗日插三种插值方法：拉格朗日插 值，分段线性值，分段线性 插值，三次样条插值。插值，三次样条插值。2.2.插值的插值的 MATLAB 实现实现及插值的应用及插值的应用。4.4.数值积分的数值积分的 MATLAB 实现实现及数值积分的应用及数值积分的应用。3什么是插值什么是插值(Interpolation)？从查函数表说？从查函数表说起起查查 函函 数数 表表标准正态分布函数表标准正态分布函数表求求 (1.114)(1.114)=0.8665 (0.8686 0.8665)0.4=0.8673插值插值插值在图像处理插值在图像处理/数控加工数控加工/外观设计等领域有重要应用外观设计等领域有重要应用4插值的基本原理插值的基本原理插值问题的提法插值问题的提法已知已知 n+1n+1个节点个节点其中其中互不相同，不妨设互不相同，不妨设求任一插值点求任一插值点处的插值处的插值节点可视为由节点可视为由产生产生,g表达式复杂表达式复杂,甚至无表达式甚至无表达式5求解插值问题的基本思路求解插值问题的基本思路构造一个构造一个(相对简单的相对简单的)函数函数通过全部节点通过全部节点,即即再用再用计算插值，即计算插值，即插值的插值的基本原理基本原理61.1.拉格朗日拉格朗日(Lagrange)(Lagrange)多项式插值多项式插值1.0 1.0 插值多项式插值多项式求三种插值三种插值方法方法71.1 1.1 拉格朗日插值多项式拉格朗日插值多项式又又（2）有唯一解，故有唯一解，故（3）与与（1）相同。相同。基函数基函数三种插值三种插值方法方法8三种插值三种插值方法方法1.2 1.2 误差估计误差估计91.3 1.3 拉格朗日插值多项式的振荡拉格朗日插值多项式的振荡Runge现象现象取n=2,4,6,8,10,计算Ln(x),画出图形三种插值三种插值方法方法Runge.m102.2.分段线性插值分段线性插值xjxj-1xj+1x0 xn计算量与计算量与n n无关无关;n n越大，误差越小越大，误差越小.三种插值三种插值方法方法11机翼下轮廓线3.3.三次样条插值三次样条插值样条函数的由来样条函数的由来飞机、船体、汽车外形等的放样（设计）飞机、船体、汽车外形等的放样（设计）细木条：样条细木条：样条123.3.三次样条插值三次样条插值数学样条（数学样条（spline）3）3）2），3）共 4n-2个方程三种插值三种插值方法方法13三次样条插值确定三次样条插值确定4 4n n个系数需增加个系数需增加 2 2个条件个条件思考1）自然边界条件的几何意义是什么？2）样条插值为什么普遍用3次多项式，而不是2或4次？三次样条三次样条插值插值14三种插值方法小结三种插值方法小结 拉格朗日插值（高次多项式插值）：拉格朗日插值（高次多项式插值）：曲线光滑；误差估计有表达式；收敛性不能保证。曲线光滑；误差估计有表达式；收敛性不能保证。用于理论分析，实际意义不大用于理论分析，实际意义不大。分段线性和三次样条插值（低次多项式插值）：分段线性和三次样条插值（低次多项式插值）：曲线不光滑（三次样条插值已大有改进）；误差估曲线不光滑（三次样条插值已大有改进）；误差估计较难（对三次样条插值）；收敛性有保证。计较难（对三次样条插值）；收敛性有保证。简单实用，应用广泛简单实用，应用广泛。其他：其他：Hermite插值、分段三次插值、二维插值等插值、分段三次插值、二维插值等根据需要，根据需要，各取所需各取所需。151.1.拉格朗日插值拉格朗日插值:自编程序自编程序,如名为如名为 lagr.m 的的M文件，文件，第一行为第一行为 function y=lagr(x0,y0,x)输入输入:节点节点x0,y0,插值点插值点x(均为均为数组，长度自定义数组，长度自定义)）；）；输出输出:插值插值y(与与x同长度数组同长度数组)）。）。应用时输入应用时输入x0,y0,x后后,运行运行 y=lagr(x0,y0,x)2.2.分段线性插值分段线性插值:已有程序已有程序 y=interp1(x0,y0,x)y=interp1(x0,y0,x,linear)3.3.三次样条插值三次样条插值:已有程序已有程序 y=interp1(x0,y0,x,spline)或或 y=spline(x0,y0,x)用MATLAB作插值计算注：注：MATLAB有样条工具箱（有样条工具箱（Spline Toolbox）16用MATLAB作插值计算为例，作三种插值的比较为例，作三种插值的比较以以 0 1.0000 1.0000 1.0000 1.0000 0.5000 0.8000 0.8434 0.7500 0.8205 1.0000 0.5000 0.5000 0.5000 0.5000 1.5000 0.3077 0.2353 0.3500 0.2973 2.0000 0.2000 0.2000 0.2000 0.2000 2.5000 0.1379 0.2538 0.1500 0.1401 3.0000 0.1000 0.1000 0.1000 0.1000 3.5000 0.0755 -0.2262 0.0794 0.0745 4.0000 0.0588 0.0588 0.0588 0.0588 4.5000 0.0471 1.5787 0.0486 0.0484 5.0000 0.0385 0.0385 0.0385 0.0385 x y y1 y2 y3 用用n=11个节个节点，点，m=21个个插值点，三插值点，三种方法作插种方法作插值，画图。值，画图。chazhi117插值的应用加工时需要加工时需要x每每改变改变0.05时的时的y值值chazhi2图1 零件的轮廓线 （x间隔0.2）表1 x间隔0.2的加工坐标x,y（图1右半部的数据）数控机床加工零件数控机床加工零件 0.0,5.00 0.2,4.710.4,4.31 0.6,3.68 0.8,3.051.0,2.50 1.2,2.051.4,1.69 1.6,1.40 1.8,1.182.0,1.00 2.2,0.862.4,0.74 2.6,0.64 模型模型 将图1逆时针方向转90度，轮廓线上下对称，只需对上半部计算一个函数在插值点的值。图2 逆时针方向转90度的结果令v=x,u=-y 18为什么要作数值积分为什么要作数值积分 许多函数许多函数“积不出来积不出来”,只能用数值方法，如只能用数值方法，如 积分是重要的数学工具，是微分方程、概率积分是重要的数学工具，是微分方程、概率论等的基础；在实际问题中有直接应用。论等的基础；在实际问题中有直接应用。对于用离散数据或者图形表示的函数对于用离散数据或者图形表示的函数，(可以先做插值然后积分；或者直接利用点做数值积分)计算积分只有求助于数值方法。计算积分只有求助于数值方法。数值数值积分积分19数数 值值 积积 分分 的的 基基 本本 思思 路路回回 忆忆 定定 积积 分分 的的 定定 义义各种数值积分方法研究的是各种数值积分方法研究的是如何取值，区间如何取值，区间如何划分，如何划分，使得既能保证一定精度，计算量又小。使得既能保证一定精度，计算量又小。n n充分大时充分大时I In n就是就是I I的数值积分的数值积分（计算功效：算得准，算得快）（计算功效：算得准，算得快）201.1.从矩形公式到梯形公式数值积分数值积分yy=f(x)xbao平均，得到梯形公式xk+1xkxk-1fk212.2.辛普森辛普森(Simpson)(Simpson)公式公式（抛物线公式）（抛物线公式）梯形公式相当于用分段线性插值函数分段线性插值函数代替每段要用相邻两小区间两小区间端点的三个函数值端点的三个函数值抛物线抛物线公式公式提高精度提高精度分段二次插值函数分段二次插值函数数值积分数值积分yy=f(x)xbaox2kf2kx2k+1x2k+2f2k+1f2k+2区间数必须为偶数区间数必须为偶数22 对对k求和求和(共共m段段)，得（复合），得（复合）辛普森公式辛普森公式：二次插值函数sk(x)2.2.辛普森辛普森(Simpson)(Simpson)公式（抛物线公式）公式（抛物线公式）23梯形公式在每小段上是用梯形公式在每小段上是用线性插值函数线性插值函数T T(x)代替代替 f(x)梯形公式梯形公式的误差估计的误差估计因为：因为：(x-xk)(x-xk+1)在在(xk,xk+1)不变号，所以：不变号，所以：24梯形公式梯形公式Tn的的误差是误差是h2阶的阶的估计估计 梯形公式梯形公式的误差的误差25同理可得：同理可得：其中其中辛普森公式辛普森公式Sn的误差是的误差是h4阶的阶的。辛普森公式的误差估计辛普森公式的误差估计26梯形公式和辛普森公式的收敛性若对若对I某个数值积分某个数值积分In有有（非零常数）非零常数）则称则称 In是是 p 阶收敛的阶收敛的。梯形公式梯形公式 2 2 阶收敛，辛普森公式阶收敛，辛普森公式 4 4 阶收敛。阶收敛。c=0:至少至少p阶收敛（超阶收敛（超p阶收敛）阶收敛）27积分步长的自动选取积分步长的自动选取选定数值积分公式后，如何确定步长选定数值积分公式后，如何确定步长h以满足给定的误差以满足给定的误差 梯形公式用二分法只要用二分法只要其中其中fk+1/2是原是原分点分点xk,xk+1的中点的中点（记记xk+1/2）的函数值的函数值且且T T2n2n可在可在T Tn n基础上计算基础上计算28高斯高斯(Gauss)(Gauss)求积公式求积公式矩形公式矩形公式(1)、(2)梯形公式梯形公式(3)辛普森公式辛普森公式(4)A Ak k是与是与f f无关的常数无关的常数代数代数精度精度设设用用(7)计算计算若对于若对于都有都有而当而当则则称称In的代数精度为的代数精度为m.Newton-Cotes方法方法29梯形公式的代数精度（考察梯形公式的代数精度（考察T1）k=1f(x)=xk=2f(x)=x2梯形公式的代数精度为梯形公式的代数精度为1辛普森公式的代数精度为辛普森公式的代数精度为330高斯公式的思路高斯公式的思路取消对节点的限制，按照代数精度最大取消对节点的限制，按照代数精度最大的原则，同时确定节点的原则，同时确定节点xk和系数和系数Ak构造求积公式构造求积公式对于对于使使G G2 2的代数精度为的代数精度为3 3确定确定31将将f(x)f(x)代入计算得代入计算得用用n个节点，个节点，Gn的代数精度可达的代数精度可达2n-1,但是需解但是需解复杂的非线性方程组，实用价值不大。复杂的非线性方程组，实用价值不大。32常 用 的 高 斯 公 式将将(a,ba,b)分小，把小区间变换为分小，把小区间变换为(-1(-1，1),1),再用再用G G2 2代数精度为代数精度为3节点加密时，原计算信息无法利用节点加密时，原计算信息无法利用33思路思路：将积分区间分小，在小区间上用：将积分区间分小，在小区间上用n不太不太 大大的的 。而在节点加密一倍时能够利用原节点的函。而在节点加密一倍时能够利用原节点的函数值，可以把区间的端点作为固定节点。数值，可以把区间的端点作为固定节点。改进的高斯公式Gauss-Lobatto求积公式求积公式 其中a,b为小区间的端点，为2n-2个参数，代数精度可达到代数精度可达到2n-3注意：实际计算中一般采用自适应方法确定步长注意：实际计算中一般采用自适应方法确定步长34用用MATLAB 作数值积分作数值积分矩形矩形公式公式Sum(x)输入数组x(即fk),输出x的和(数)cumsum(x)输入数组x,输出x的依次累加和(数组)梯形梯形公式公式trapz(x)输入数组x,输出按梯形公式x的积分(单位步长)trapz(x,y)输入同长度数组 x,y,输出按梯形公式y对x的积分(步长不一定相等)35用用MATLAB 作数值积分作数值积分辛普森公式辛普森公式quad(fun,a,b,tol,trace)I,fn=quad()用自适应辛普森公式计算tol为绝对误差，缺省时为10-6Gauss-Lobatto公式公式quadl(fun,a,b,tol,trace)I,fn=quadl()用自适应Gauss-Lobatto公式计算 tol为绝对误差，缺省时为10-6注意：注意：fun.m中应以自变量为矩阵的形式输入中应以自变量为矩阵的形式输入(点运算点运算)36矩形域上计算二重积分的命令：矩形域上计算二重积分的命令：dblquad(fun,xmin,xmax,ymin,ymax,tol)广义积分、二重和三重积分长方体上计算三重积分的命令：长方体上计算三重积分的命令：triplequad(fun,xmin,xmax,ymin,ymax,zmin,zmax,tol)注：注：fun是被积函数，本身可以有自己的参数是被积函数，本身可以有自己的参数广义积分：广义积分：通过分析和控制误差，转换成普通积分通过分析和控制误差，转换成普通积分quadv(fun,a,b,tol,trace)向量值积分：向量值积分：37用用MATLAB 作数值积分作数值积分例例.计算计算1 1）矩形公式和梯形公式）矩形公式和梯形公式将将(0,(0,/4)100/4)100等分等分2 2）辛普森公式和辛普森公式和Gauss-Lobatto公式公式精确、方便精确、方便无法计算用数值给出的函数的积分无法计算用数值给出的函数的积分Jifen1a.mJifen1a.mJifen1b.mJifen1b.m精确值为精确值为38数值积分的应用数值积分的应用实例实例人造卫星轨道长度人造卫星轨道长度轨道长度轨道长度yxo 近地点s1=439km,远地点s2=2384kms1s2地球半径r=6371kmr需要作数值积分需要作数值积分39s1=439km,s2=2384km,r=6371kmyxos1s2rs1s2yxoracb数值积分实例数值积分实例 人造卫星轨道长度人造卫星轨道长度40用梯形公式和辛普森公式计算用梯形公式和辛普森公式计算只将区间只将区间5等分，梯形公式就给出很好的结果等分，梯形公式就给出很好的结果轨道长度轨道长度 L L=4.8707 104千米千米数值积分实例数值积分实例 人造卫星轨道长度人造卫星轨道长度Jifen2.mJifen2.m41布置实验布置实验目的目的1、掌握用、掌握用MATLAB计算拉格朗日、分段线性、三次计算拉格朗日、分段线性、三次样条三样条三 种插值的方法，改变节点的数目，对三种种插值的方法，改变节点的数目，对三种插值结果进行初步分析。插值结果进行初步分析。2、掌握用掌握用MATLAB及梯形公式、辛普森公式计算及梯形公式、辛普森公式计算数值积分。数值积分。3、通过实例学习用插值和数值积分解决实际问题。通过实例学习用插值和数值积分解决实际问题。内容内容实验练习实验练习 1010；1212；42