数值微分和数值积分

2015/1/5数值计算方法1数值微分和数值积分浙江大学控制系微分和积分(DifferentiateandIntegrate)2015/1/5数值计算方法20()()()()lim()iiiixbafxxfxyxxfxxfxdydxxIfxdx微分和积分Newton-Leibniz公式2015/1/5数值计算方法3)()()(aFbFdxxfba微分和积分函数类型简单的连续函数，如多项式、指数函数或三角函数很难或完全无法直接微分或积分的复杂函数列表型函数，仅给出了一系列离散点及相应的函数值，如实验或现场数据利用插值或拟合得到的函数2015/1/5数值计算方法4等面积微分2015/1/5数值计算方法5网格逼近积分2015/1/5数值计算方法6f(x)abx积分中值定理2015/1/5数值计算方法7))(()(abfdxxfbaf(x)abx数值微分与积分数值微分有限差商逼近导数如果数据存在误差，使用曲线拟合技术构造光滑拟合曲线，对拟合曲线微分数值积分积分中值定理用区间[a,b]上一些点xk处的“高度”f(xk)的加权平均值，作为平均高度f()的近似值2015/1/5数值计算方法8()()()bafxdxfba0()()nbkkakIfxdxAfx求积系数Ak=k(ba)，k是求积节点权系数。数值微分与积分对于复杂连续函数由函数表达式生成离散的数据列表在离散点的基础上，利用数值方法计算数值微分和积分2015/1/5数值计算方法9数值积分与微分的稳定性一般数值积分过程是稳定的，所得的解精确度也较高。数值微分时，因近似多项式的曲线斜率可能和给定函数f(x)的曲线斜率有很大不同，特别是当f(x)在给定区间内变化比较大时更是这样。这样就使得数值微分的解不稳定，并且精度也较差。2015/1/5数值计算方法10f(x)x0abP(x)积分过程对数据点进行求和，正的和负的随机误差倾向于相互抵消；相反，微分过程是相减的，正的和负的随机误差倾向于相加。本章内容数值微分差商近似插值型数值微分数值积分Newton-Cotes积分龙贝格（Romberg）积分高斯（Gauss）求积公式2015/1/5数值计算方法11有限差商近似向前差商向后差商中心差商2015/1/5数值计算方法12000()()()fxhfxfxh000()()()()()()2!fxhfxhRxfxfOhh000()()()fxfxhfxh000()()()()()()2!fxfxhhRxfxfOhh000()()()2fxhfxhfxh00022212()()()()2[()()]()()126fxhfxhRxfxhhhfffOh有限差商近似的误差h越小，误差越小，但同时舍入误差增大最佳步长确定的事后估计法设D(h)，D(h/2)分别为步长为h，h/2的差商公式2015/1/5数值计算方法13()()()()(/2)(/2)fxDhOhfxDhOh()()()2()(/2)(/2)fxDhOhfxDhOh()()2()2(/2)fxDhfxDh()(/2)()(/2)fxDhDhDh)2()(hDhD时的步长h/2就是合适的步长增加泰勒级数展开式中的项数可以提高精度高精度微分公式2015/1/5数值计算方法14221()4()3()()()2iiiifxfxfxfxOhh21()()()()2iiiifxfxfxfxhh21()()()()()2iiiifxfxfxfxhOhh212()2()()()()iiiifxfxfxfxOhh21212()()()2()()()()2iiiiiifxfxfxfxfxfxhOhhh增加二阶导数项可以将精度提高到O(h2)数值微分——例例：利用有限差商估计函数在x=0.5处的导数值（真值为f’(0.5)=-0.9125）解：h=0.5，前向差商后向差商中心差商h=0.252015/1/5数值计算方法15432()0.10.150.50.251.2fxxxxxixif(xi)i-101.2i0.50.925i+11.00.20.20.925(0.5)1.4558.9%0.5tf0.9251.2(0.5)0.5539.7%0.5tf0.21.2(0.5)1.09.6%1tfixif(xi)i-10.251.10351563i0.50.925i+10.750.63632813(0.5)1.15526.5%tf(0.5)0.71421.7%tf(0.5)0.9342.4%tf数值微分——例取h=0.25，利用高精度微分公式前向差商后向差商中心差商2015/1/5数值计算方法16ixif(xi)i-201.2i-10.251.10351563i0.50.925i+10.750.63632813i+21.00.2(0.5)0.8593755.82%tf(0.5)0.8781253.77%tf(0.5)0.91250%tf本章内容数值微分差商近似插值型数值微分数值积分Newton-Cotes积分龙贝格积分高斯求积公式2015/1/5数值计算方法17插值型数值微分用插值函数的导数近似为原函数的导数2015/1/5数值计算方法18()()()()kknfxPx(1)1()()()(1)!nnnfRxxn(1)()1()()()(1)!knknnkdfRxxdxn,0,1,2,,jjxfxjn)()()(xRxLxfnn)()()(xRxNxfnn插值型数值微分k=1未知，当xxi时，无法利用上式估计误差，若求某个插值节点的导数值，有特例：n=1，插值节点为x0和x1，若记h=x1-x0——一阶微分两点公式2015/1/5数值计算方法19(1)(1)11()()()()()(1)!(1)!nnnnnxfdRxxfnndx(1)1()()()()(1)!nininiffxPxxn100()()'()()2fxfxhfxfh101()()'()()2fxfxhfxfh在x0处的向前差商在x1处的向后差商插值型数值微分特例：n=2，插值节点为xk=x0+kh(k=0,1,2)——一阶微分三点公式2015/1/5数值计算方法202(3)00122(3)1022(3)20121()3()4()())()231()()()()261()()4()3()()23hfxfxfxfxfhhfxfxfxfhhfxfxfxfxfh中心差商三点向前差商三点向后差商基于泰勒级数的公式等价于通过对数据点进行插值求微分数值微分用三次样条插值的导数近似被插值函数的导数，效果相当好在实际问题中使用哪种公式要视具体问题而定。有时三种公式都要用到，给定一个列表函数，对函数中间各点都可使用精度较高的中心差分公式，但起始点只能使用前差公式，终点则使用后差公式。一般情况下，三点公式比两点公式准确，步长越小结果越准确。但当余项中的高阶导数无界或计算过程中的舍入误差超过截断误差时，这个结论不成立。2015/1/5数值计算方法21本章内容数值微分差商近似插值型数值微分数值积分Newton-Cotes积分龙贝格（Romberg）积分高斯（Gauss）求积公式2015/1/5数值计算方法22数值积分定义数值积分是离散点上的函数值的线性组合两个问题：求积系数如何选取节点可以自由选取，取什么点好2015/1/5数值计算方法230()()()nbiinaiIfxdxAfxIf()()nnIIfRf称为求积系数，与f(x)无关，与积分区间和积分点有关求积公式的余项称为求积节点，与f(x)无关代数精度如果一个求积公式对任何次数不超过m次的代数多项式都准确成立，而对于m+1次的代数多项式不能准确成立，即R(f)0，则称求积公式In(f)的代数精度是m梯形公式的代数精度m=1对不高于1次的代数多项式都准确成立对2次以上的代数多项式存在误差一个求积公式的代数精度越高，就能对更多的被积函数准确或较准确地成立2015/1/5数值计算方法24插值型求积公式用插值函数的积分，作为数值积分代数精度2015/1/5数值计算方法25)()()()()()(00inibaiibaniibannxfdxxldxxfxldxxLfIAi(1)1()()()(1)!nnnfRxxn(1)1()()()()()(1)!nbbnnnaafIfIfRxdxxdxn(1)()0,(),nkfxfxxkn至少n阶代数精度本章内容数值微分差商近似插值型数值微分数值积分Newton-Cotes积分龙贝格（Romberg）积分高斯（Gauss）求积公式2015/1/5数值计算方法26Newton-Cotes积分2015/1/5数值计算方法27Newton-Cotes积分2015/1/5数值计算方法28Newton-Cotes积分采用拉格朗日插值多项式Pn(x)来逼近f(x)，将积分区间n等分，取分点为求积节点，并作变量替换2015/1/5数值计算方法29niihaxnabhi,,0,,,xath00(1)(1)(1)()()!()!(1)(1)(1)(1)(1)()!()!bniinianintttititnAlxdxhdtininhtttititndtininb-a与步长h无关，可以预先求出——Cotes系数)(niC()()niiAbaC梯形公式n=12015/1/5数值计算方法302121)1(10)1(110)1(0dttCdttC111()()()()()22()()()2Ifbafabafbfafbba13()()()()2!()()()2!()()12babafEfxaxbdxfxaxbdxbafI=(b-a)平均高度梯形公式——例例：利用梯形公式计算函数在[00.8]上的积分（精确值为1.640533）解：f(0)=0.2f(0.8)=0.232估计误差：用被积函数二阶导数的积分平均值代替f”()2015/1/5数值计算方法312345()0.225200675900400fxxxxxx00.10.20.30.40.50.60.70.800.511.522.533.5411(0.80)(0.20.232)0.17282I89.5%t0.8230(4004050108008000)()600.80xxxdxfx1.467733tE31(60)(0.8)2.5612aESimpson公式2015/1/5数值计算方法32Simpson公式n=22015/1/5数值计算方法3361)1(4164)2(2161)2)(1(4120)2(220)2(120)2(0dtttCdtttCdtttC2141()()()()()()()6626()4()()32baIfbafabafbafbhbafaffbSimpson1/3法则I