數(shù)值微分的計算方法

1、數(shù)值微分的計算方法內(nèi)容摘要 求解數(shù)值微分問題，就是通過測量函數(shù)在一些離散點上的值，求得函數(shù)的近似導(dǎo)數(shù)。本文就所學(xué)知識，歸納性地介紹了幾種常用的數(shù)值微分計算方法。并舉例說明計算，實驗結(jié)果表明了方法的有效性。關(guān)鍵詞數(shù)值微分 Taylor 展開式Lagrange 插值三對角矩陣引言：數(shù)值微分即根據(jù)函數(shù)在一些離散點的函數(shù)值，推算它在某點的導(dǎo)數(shù)或高階導(dǎo)數(shù)的近似值的方法。常見的可以用一個能夠近似代替該函數(shù)的較簡單的可微函數(shù)(如多項式或樣條函數(shù)等)的相應(yīng)導(dǎo)數(shù)作為能求導(dǎo)數(shù)的近似值，由此也可導(dǎo)出多點數(shù) 值微分計算公式。當(dāng)函數(shù)可微性不太好時，利用樣條插值進(jìn)行數(shù)值微分要比多項式插 值更適宜。I. Taylor 展開

2、式方法理論基礎(chǔ)：Taylor展開式2nf(X)=f(Xo 廣(X冷/ J )+2f"(X。)*+fCX°我們借助Taylor展開式，可以構(gòu)造函數(shù) f (x)在點x= X。的一階導(dǎo)數(shù)和二階導(dǎo)數(shù)的數(shù)值微分公式。取步長 h。則(1)(2)(3)(4)f (Xo h) = f (x。) hf'(x。) h f ( i) i (X。，X。 h) 2所以f (Xo h)f (Xo) h .f (X。)- f ( i) i(號。h)h2同理2 (Xo h,X。)hf (Xo -h) =f (X。)-hf (X。)萬 f ( 2)f'(5。)-(尚-h) §f&#

3、39;'(偵 2 (X。一h,X。)h2式 和式(4)是計算f'(X。)的數(shù)值微分公式，其截斷誤差為0(h),為提高精度，將Taylor展開式多寫幾項f(xo h) = f(Xo) hf2 3 4(X。)3 f (X。)X f (X。)板 f ( 1)2624, h2f (X。-h) = f (X。)一hf (X。)23 4f (X。)- ： f (X。)h f(4)( 2)6241 (乂。，*。 h)2 (X。-h, X。)兩式相減得-2(5)f ( x h) - f(Xnh) h'"4f (X。)= f(X° h) f(X° h)_ L

4、 f (xo) 0(h )2h62上式為計算f (X。)的微分公式，其截斷誤差為 0(h),比式(2)和(4)精度寄兩式相加，如果f (4)(X)ECX。一 h,x。+h,則有2(6)f''(x)_f(x)h)_2f(X0)f(X0_h) h f(4)( f(Xo) 一h2一12 f J)E(xoh,x°+h)2式(6)是計算f (Xo)的數(shù)值微分公式，其截斷誤差為O(h2)。例1 .設(shè)函數(shù)f(x)= lnx,X0= 2,h= 0.1，試用數(shù)值微分公式計算f (2)的值。解 由式(2)、式和式(5)分別計算結(jié)果為ln 2.1ln 2f (2) 0.48790.1ln

5、2 - ln 1.9 - 一f (2)：、: 0.51290.1ln2.1_ln1.9f (2) 0.50040.2與真值f'(2)= 0.5相比，式(5)計算的結(jié)果精度較高。2. 數(shù)值微分的Lagrange插值方法設(shè)函數(shù)y= f (x)具有n+ 1個實驗數(shù)據(jù)：(為，f (xi)(i = 0,1,2，n),我們希望估討'，、f (x)的值，特別x =x時，估計f (x)的值?；诓逯捣椒ǖ臄?shù)值微分做法是，由已知(x, f (x)(i = 0,1，,n)建立Lagrange插值多項式或Newton插值多項式(這里以 Lagrange插值方法為例)，即f (x)= Ln (x)+

6、尺(x)于是f x = Ln xRn xx= x時，有f (xi) =Ln(xi)Rn(xi)(i =0,1,2,., n)其中'f (n 1)( )，Rn(x) = n1(xi)(n 1)!略去誤差項有f (x) ： Ln(xi)實際運用中，等距節(jié)點更為常見。設(shè), b- ah =，為=na + ih (i = 0,1,2，,n), x = a + th ,于是有f (x) = y° 十；偵 +")M 十.EaM + Rn(X)所以d , tf (為)=?y0y0 -.-1!()n(1 -j)(n 1)!j=0jTdtn (n 1) h fn+ Ktf.Jt n 1

7、) ,：ny0tgn!3. 多點數(shù)值微分公式由于高階插值的不穩(wěn)定性, 型求導(dǎo)公式。(1)兩點公式(n=1)實際應(yīng)用時多采用n=1,2,4的兩點、三點和五點等多點插值f'(X0) =：(yi -y0)h、1 ,、f (X1) =(y -y°) h_ ,h _ R (X0) 2 f ()R(X1)業(yè) f''()2(2) 三點公式(n=2),',、1 ，、f (冷)=二(3y° 4乂 y?)2h一，,、1 .f (X1) = = (y。+y2)2h,',1f (X2) =(y° 4山 +3y2)2hR(X0)=§f(3)

8、()R(X1) =-項 f(3)()62R(X0) = §f(3)()3(8)五點公式(n=4)r,',、1,f (x°)(-25y° 48y -36y2 16y3 -3y4)12h.',、1,f (X1)(-3y0 -10y1 18y2 -6y3 - y4)12h，、1,f (X2)(y° -8y1 8y3 -y4)12h,'1f (X3)(T0 6貝18y2 10y3 3y4)12h，、1,f (X4)(3y° -16y" 36y48y3 25y4)12hR(X0)=gf(5)()5R(xJA* )20R(

9、X2h：f(5)() 30R(X3)=-£f(5)()20R(X°)=£f(5)()5(9)例2.設(shè)f (x) = ln x，取h = 005 ,分別用三點公式和五點公式計算f' (2)的近似值。解：由式(8)有1f (2). .(-3f (2) 4 f (2.05) f (2.10) =0.4998028612 0.05f'(2) =1(f (1.95) - f (2.05) =0.5001042052 0.05 1f (2). . (f (1.90) -4f (1.95) - 3f (2) =0.4997793762 0.05由式(9)有1f

10、(2) =1(f (1.90) _8 f (1.95) 6f (2.05) _f (2.10) =0.49999984312 0.05t . _f. t . 一' 一、_. .與真值f (2)= 0.5相比，三點公式已有相當(dāng)滿意精度，而五點公式的結(jié)果是十分滿意的。4. 數(shù)值微分的隱式格式前述的Taylor 展式法和數(shù)值微分格式均稱為顯式格式，即直接由已知的 ， 、 . . '， 、f (Xi)(i = 0,1,2，n)，經(jīng)過適當(dāng)?shù)乃阈g(shù)四則運算，立即可得 f (X)的近似值。顯式格式優(yōu)點是計算方便，工作量小，缺點是數(shù)值不穩(wěn)定。為克服后一缺點，隱式格式常常具有數(shù)值穩(wěn)定性。數(shù)值微分的

11、隱式格式建立方法常用通過Taylor展開式方法或數(shù)值積分方法等不同途徑。首先我們用Taylor展開式方法來推導(dǎo)數(shù)值微分的隱式格式。由式和式(6)，我們用&代替X0得：1h'”一 4f(XjF皿 h)-f(Xk 頂-百 f ("0(h)f''(Xk) =f(Xk h) -2f (Xk) f(Xk h) f (4)()h212所以也有'''1'''h2(5)-f同二仍心為-小口叢叫-"？()將最后f (xk)表達(dá)式代入f (Xk)表達(dá)式可得114f (Xk)Fg 枷-皿-h)飛仃(Xk +h) -2

12、f (f (Xk -h) +o(h)略去誤差項o(h4)，并且mk表示f (Xk)的近似值，則有,3 、mk 1 - 4mk mk 1 f (Xk 1) - f (Xk 盤(k =1,2,., n)h同樣，用數(shù)值積分方法也可推導(dǎo)數(shù)值微分的隱式格式。Xk 1Xk ，f(X)dX=f'(Xk )dX(k = 1,2,., n)'、X _Lbb -aa b.根據(jù) simpson 公式f(X)dX =f (a) +4f( 一)+ f (b)有 a 、，622hf (Xk 1) - f (Xk)=;mk4m" mk 1, (k = 1,2,., n)6,3，mk4mk mk 1

13、 =f g 1) - f (Xj)(k =1,2,., n)h上式是關(guān)于n+ 1個未知量m0,m,mn的n- 1個方程，如果已知 m0 =3 f (X°),mn = f (Xn)，記 dk = f (Xf) f (Xj),故有萬程組：h4 1一 m1-d m°1 4 1m2d2141m3d3+ + +a， > ， 14 1mudn:. 1 4 一：mn一idn-mn_(10)這就是求m1,典，mn-1的線性方程組。由于系數(shù)矩陣是嚴(yán)格對角占優(yōu)的三對角矩陣，因此非奇異，解存在且唯一，可由追趕法求解，且數(shù)值穩(wěn)定。例3.設(shè)函數(shù)f (x)= lnx在節(jié)點x= 1.5, 1.6,

14、 1.7, 1.8, 1.9上的函數(shù)值，及-'. f (1.5), f (2.0)如表1所示。試用數(shù)值微分隱式格式式(10)和式(11)求出相應(yīng)節(jié)點上一階導(dǎo)數(shù)的近似值。表1例3節(jié)點數(shù)據(jù)xf(x)f ' (x)1.50.4054651080.6666666671.60.4700036291.70.5306282511.80.5877866641.90.6418538862.00.6931471820.500000000解：由式(10)有面m21 m3七郵一3.75489429 -0.6666666673.533491053.33676905316081554 -0.5000000

15、00-解方程組得:m =0.62499828611483m2 =0.58823447854067ms =0.55555484972249m4 =0.52631517256938.1 與f (x)=在各相應(yīng)節(jié)點上數(shù)值相比，上述結(jié)果約有五位有效數(shù)字，精度較高。h越小,精度越高。附追趕法程序：%追趕法%定義三對角矩陣 A的各組成單元。方程為 Ax=d% a為主對角線下面的次對角線，b為主對角線，c為主對角線上面的次對角線，d為等式右端向量。% A=4100%1410%0141%0014function x=zhuiganfa(a,b,c,d)a=0 1 1 1;c=1 1 1;b=4 4 4 4;d

16、=3.75489429-0.666666667,3.53349105,3.33676905,3.16081554-0.5;n=length(b);u0=0;y0=0;a(1)=0;L(1)=b(1)-a(1)*u0;y(1)=(d(1)-y0*a(1)/L(1);u(1)=c(1)/L(1);for i=2:(n-1)L(i)=b(i)-a(i)*u(i-1);y(i)=(d(i)-y(i-1)*a(i)/L(i);u(i)=c(i)/L(i);endL(n)=b(n)-a(n)*u(n-1);y(n)=(d(n)-y(n-1)*a(n)/L(n);x(n)=y(n);for i=(n-1):-1:1x(i)=y(i)-u(i)*x(i+1);end運行結(jié)果：ans =0.624998286114830.588234478540670.555554849722490.52631517256938>> x=1.5:0.1:2.0;0.588234478540670.55555484972249y=log(x);dy=1/1.50.624998286114830.52631517256938 1/2.0;plot(x,y,'b-'