西安交通大学：《计算物理学》研究生课程教学资源（课件讲稿）第一章 基本数学运算（主讲教师：蒋臣威）

本章内容XI'ANJIAOTONGUNIVERSITY数值微分2数值积分3方程求根分子振动的半经典量子化交道大学

XI’AN JIAOTONG UNIVERSITY 本章内容 1 数值微分 2 数值积分 3 方程求根 4 分子振动的半经典量子化5

1.数值微分XI'ANJIAOTONGUNIVERSITY为什么要学习数值微分？我们碰到的函数往往没有解析形式，例如可能是如下数表形式如利用计算所得势能曲线上的点求梯度力。0.50.10.20.30.4x6844.5f(x)8.5这就需要借助于数值微分。当函数f(x）过于复杂时，也可借助数值微分。更重要的，数值微分是其它很多数值方法的基础

XI’AN JIAOTONG UNIVERSITY 1. 数值微分 为什么要学习数值微分？ 我们碰到的函数往往没有解析形式，例如可能是如下数表形式 如利用计算所得势能曲线上的点求梯度力。 x 0.1 0.2 0.3 0.4 0.5 f(x) 4 4.5 6 8 8.5 这就需要借助于数值微分。 当函数f(x)过于复杂时，也可借助数值微分。 更重要的，数值微分是其它很多数值方法的基础

XIANJIAOTONGUNIVERSITY微积分中，关于导数的定义如下：f(x+h)- f(x-h)f(x+ h)- f(x) - lim f(x)- f(x-h)f(x)= lim Tinlim2hhhh-→0h→0h→0自然，而又简单的方法就是，取极限的近似值，即差商（差分）向前差分f(xo +h)- f(x)f(xo)~)h由Taylor展开h?f(xo+h)= f(x)+hf(x-2!h?"(E),xo ≤=≤x +h= f(xo)+ hf '(xo)2!h= f(x0)= f(xo +h)-f(x0)f"(), x ≤≤ x +hh21

XI’AN JIAOTONG UNIVERSITY 微积分中，关于导数的定义如下： 000 ( ) () () ( ) ( ) ( ) '( ) lim lim lim hhh 2 fx h fx fx fx h fx h fx h f x ®®® hh h + - - - + - - === 自然，而又简单的方法就是，取极限的近似值，即差商（差分） 0 0 0 ( ) () '( ) fx h fx f x h + - » 由Taylor展开 2 0 00 0 2 0 0 00 ( ) () () () 2! ( ) '( ) ''( ), 2! h f x h f x hf x f x h f x hf x f x x h x x += + + + ¢ ¢¢ ××× =+ + £ £ + 向前差分 0 0 0 00 ( ) () '( ) ''( ), 2! fx h fx h fx f x x h h x x + - Þ = - ££ +

XI'ANJIAOTONGUNIVERSITY因此，有误差hf(xo +h)- f(xo)R(x)= f(xo)- ()= 0(h)h2!向后差分f(xo)- f(xo -h)f'xh由Taylor展开h?f(xo -h)= f(xo)- hf'(xo+2!h?F"(E),x≤=≤x。 + h= f(xo)- hf'(xo2!爱通大学nf(xo)- f(x。-h)"(5),xo ≤5≤x。 + h=f'(x2h

XI’AN JIAOTONG UNIVERSITY 因此，有误差 0 0 0 ( ) () ( ) '( ) ''( ) ( ) 2 ! fx h fx h Rx f x f O h x h + - = - = - = 0 0 0 () ( ) '( ) fx fx h f x h - - » 2 0 00 0 2 0 0 00 ( ) () () () 2! ( ) '( ) ''( ), 2! h f x h f x hf x f x h f x hf x f x x h x x - = - ¢ ¢¢ + -××× = - + £ £ + 向后差分 由Taylor展开 0 0 0 00 () ( ) '( ) ''( ), 2 fx fx h h fx f x x h h x x - - Þ = + £ £ +

XI'ANJIAOTONGUNIVERSITY因此，有误差f(xo)-f(x-h) hR(x) = f'(xof "() = 0(h)2!h中心差分f(xo + h)- f(xo - h)f'(x.2 h由Taylor展开hh?f(x +h)= f(x)+hf'(x3!2!h3h?(5),x ≤Ei ≤ x + h= f(xo)+ hf '(xo3!2!h3h2f(xo - h)= f(xo)- hf'(xo+：3!2!复通大学h2h3= f(xo)- hf'(x(52),x0 - h≤52 ≤xo2!3!

XI’AN JIAOTONG UNIVERSITY 0 0 0 () ( ) ( ) '( ) ''( ) ( ) 2 ! fx fx h h Rx f f O x h h x - - = - = = 0 0 0 ( )( ) '( ) 2 fx h fx h f x h + - - 中心差分 » 因此，有误差 2 3 0 00 0 0 2 3 0 0 0 10 1 0 ( ) ( ) ( ) ( ) '''( ) 2! 3! ( ) '( ) ''( ) '''( ), 2! 3! h h f x h f x hf x f x f x h h f x hf x f x f x x h x x += + + + + ¢ ¢¢ ××× =+ + + £ £ + 由Taylor展开 2 3 0 00 0 0 2 3 0 0 0 20 2 0 ( ) ( ) ( ) ( ) '''( ) 2! 3! ( ) '( ) ''( ) '''( ), 2! 3! h h f x h f x hf x f x f x h h f x hf x f x f x h x x x - = - ¢ ¢¢ + - + ××× = - + - -£ £

XT'ANJIAOTONGUNIVERSITY上面两式相减得：h3h3f(xo + h)- f(x。- h) = 2hf'(3!3!h?J(x)= f(o +h)-f(xo -h)2.2h12f(xo +h)-f(xo-h)2h因此，有误差f(xo +h)- f(x -h)R(x) = f'(xo2hh?h2f "() = O(h?f"(S)+ f"()]=126夷道大学

XI’AN JIAOTONG UNIVERSITY 3 3 00 0 1 2 ( ) ( ) 2 ( ) '''( ) '''( ) 3! 3! h h f x h f x h hf x f f + - - = + + ¢ x x 上面两式相减得： 2 2 0 0 0 12 0 0 ( )( ) ( ) [ '''( ) '''( )] 2 12 12 ( )( ) 2 fx h fx h h h fx f f h fx h fx h h x x + - - ¢ = - + + - - » 因此，有误差 0 0 0 2 2 2 1 2 ( )( ) ( ) '( ) 2 [ '''( ) '''( )] '''( ) ( ) 12 6 fx h fx h Rx f x h h h f f fO xx x h + - - = - = - + = - =

XI'ANJIAOTONGUNIVERSITY五点公式hf"(x)+O(h*)f = f(x。±h)= f(xo)±hf'(xo)++J+2 = f(xo ±2h)= f(xo)±2hf'(xo)+2h'f"(xo)±=h'f"(xo)+O(h*)3h3分别相减f+1- f-{ = 2hf'(xo)+"- f"(xo)+O(hs)3Iα- L.= 4h(c)+号nf"(n)+(h)38(f+1 - f-,)-(f+2 - f-,)=12hf'(x)+ O(h')[8(f{-f-1)-(f+2-f-2)]+O(h") 五点公式X12h通大学虽然精度更大，但是计算量也更天！通常用三点公式已足够

XI’AN JIAOTONG UNIVERSITY 23 4 20 0 0 0 0 4 ( 2) ( ) 2 ( ) 2 ( ) ( ) ( ) 3 f f x h f x hf x h f x h f x O h ± º ±= ± + ± + ¢ ¢¢ ¢¢¢ 2 3 4 10 0 0 0 0 ( ) () () () () () 2 6 h h f f x h f x hf x f x f x O h ± º ±= ± + ± + ¢ ¢¢ ¢¢¢ 分别相减 3 5 11 0 0 2 () () () 3 h f f hf x f x O h + - - = + + ¢ ¢¢¢ 3 5 22 0 0 8 4 () () () 3 f f hf x h f x O h + - - = + + ¢ ¢¢¢ 5 11 2 2 0 8( ) ( ) 12 ( ) ( ) f f f f hf x O h + - + - - = + ¢ 4 0 11 22 1 ( ) [8( ) ( )] ( ) 1 2 f x f f f f Oh h + - + - ¢ = - + 五点公式 虽然精度更大，但是计算量也更大！通常用三点公式已足够。 五点公式

XI'ANJIAOTONGUNIVERSITY表一.不同方法计算d(sinx)/dxlx=1=0.540302的误差ForwardBackwardSymmetricSymmetrich3-point2-point2-point5-point0.500000.0010920.0222330.228254-0.1837890.200000.0035950.087461-0.0802720.0000280.100000.0008990.042938-0.0411390.0000010.050000.0002250.021258-0.0208080.0000000.020000.0000370.008453-0.0083800.0000010.010000.0000100.004224-0.0042040.0000020.005000.0000100.002108-0.0020880.0000060.00200-0.0000140.000820-0.000848-0.0000170.001000.000403-0.000431-0.000019-0.0000140.000500.000105-0.0001930.0001150.0004030.00020-0.000163-0.000014-0.000312-0.0001880.00010-0.000312-0.000312-0.000312-0.0004110.000050.0002840.0014760.000681-0.0009080.000873见Matlab程序0.000880chap1 sinx_differential.m

XI’AN JIAOTONG UNIVERSITY 表一. 不同方法计算 d(sin x)/dx |x=1 = 0.540302 的误差 见Matlab程序 chap1_sinx_differential.m

XI'ANJIAOTONGUNIVERSITY注意不同公式的精度：五点>三点>两点>注意误差随步长h的减小先减小再增大舍入误差微分公式涉及两个很接近的f值相减。当步长过小时，计算机的舍入误差会使导数计算不准确。合适步长设D(h),D(h/2)分别为步长为h,h/2的差商公式。则2h)120<8麦蔬大学时的步长h/2就是合适的步长

XI’AN JIAOTONG UNIVERSITY Ø 不同公式的精度：五点>三点>两点 Ø 注意误差随步长 h 的减小先减小再增大 注意 舍入误差 微分公式涉及两个很接近的 f 值相减。当步长过小时， 计算机的舍入误差会使导数计算不准确。 设D(h),D(h/2)分别为步长为h,h/2的差商公式。则 合适步长 () ( ) 2 h Dh D- < e 时的步长h/2就是合适的步长

高阶导数XT'ANJIAOTONGUNIVERSITY二阶导数h?h3(x) +0(h4)fi = f(xo +h)= f(xo)+hf'26h?h30)+0(h*)f-, = f(xo -h)= f(xo)-hf'(xXo26F"(x)= f(x +h)-2f(0)+ (-h)+O(h2)h?这也容易从二阶导数的定义直接得出（利用一阶导数向后差分公式）f'(xo +h)- f(xo) _ [f(x +h)-f(x)] / h-[f(x)-f(x-h)]/ hf"(x)=hhf(x +h)-2f(x)+ f(x-h)麦通大学h?更高阶的导数以此类推尝试给出一个三阶导数的公式

XI’AN JIAOTONG UNIVERSITY 高阶导数 二阶导数 更高阶的导数以此类推 这也容易从二阶导数的定义直接得出（利用一阶导数向后差分 公式） 2 3 4 10 0 0 0 0 2 3 4 10 0 0 0 0 ( ) () () () () () 2 6 ( ) () () () () () 2 6 h h f f x h f x hf x f x f x O h h h f f x h f x hf x f x f x O h - º += + + + + ¢ ¢¢ ¢¢¢ º - = - ¢ ¢¢ ¢¢¢ + - + 0 2 0 0 2 0 ( ) 2( ) ( ) () () fx h fx fx h fx O h h + - + - ¢¢ = + 0 0 0 0 00 0 0 00 2 '( ) '( ) [ ( ) ( )]/ [ ( ) ( )]/ ( ) ( ) 2( ) ( ) f x h f x fx h fx h fx fx h h f x h h fx h fx fx h h + - + - - - - ¢¢ = = + - + - = 尝试给出一个三阶导数的公式