Giáo trình C - Chương 12: Tính gần đúng đạo hàm và tích phân xác định

204
Ch−¬ng 12 : TÝnh gÇn ®óng ®¹o hµm vµ tÝch ph©n x¸c ®Þnh 
§1. §¹o hµm Romberg 
§¹o hµm theo ph−¬ng ph¸p Romberg lµ mét ph−¬ng ph¸p ngo¹i suy ®Ó x¸c ®Þnh ®¹o 
hµm víi mét ®é chÝnh x¸c cao . Ta xÐt khai triÓn Taylor cña hµm f(x) t¹i (x+h) vµ (x-h) : 
⋅⋅⋅++′′′+′′+′+=+ )x(f
!4
h
)x(f
!3
h
)x(f
2
h
)x(fh)x(f)hx(f )4(
432
 (1) 
⋅⋅⋅−+′′′−′′+′−=− )x(f
!4
h
)x(f
!3
h
)x(f
2
h
)x(fh)x(f)hx(f )4(
432
 (2) 
Trõ (1) cho (2) ta cã : 
⋅⋅⋅++′′′+′=−−+ )x(f
!5
h2
)x(f
!3
h2
)x(fh2)hx(f)hx(f )5(
53
 (3) 
Nh− vËy rót ra : 
⋅⋅⋅−−′′′−−−+=′ )x(f
!5
h
)x(f
!3
h
h2
)hx(f)hx(f
)x(f )5(
42
 (4) 
hay ta cã thÓ viÕt l¹i : 
[ ] ⋅⋅⋅++++−−+=′ 664422 hahaha)hx(f)hx(fh21)x(f (5) 
trong ®ã c¸c hÖ sè ai phô thuéc f vµ x . 
Ta ®Æt : 
)]hx(f)hx(f[
h2
1)h( −−+ϕ = (6) 
Nh− vËy tõ (5) vµ (6) ta cã : 
 ⋅⋅⋅−−−−′=ϕ= 664422 hahaha)x(f)h()1,1(D (7) 
⋅⋅⋅−−−−′=⎟⎠
⎞⎜⎝
⎛ϕ=
64
h
a
16
h
a
4
h
a)x(f
2
h
)1,2(D
6
6
4
4
2
2 (8) 
vµ tæng qu¸t víi hi = h/2
i-1 ta cã : 
⋅⋅⋅−−−−′=ϕ= 6i64i42i2i hahaha)x(f)h()1,i(D (9) 
Ta t¹o ra sai ph©n D(1,1) - 4D(2,1) vµ cã : 
⋅⋅⋅−−−′−=⎟⎠
⎞⎜⎝
⎛ϕ−ϕ 6644 ha16
15
ha
4
3
)x(f3
2
h
4)h( (10) 
Chia hai vÕ cña (10) cho -3 ta nhËn ®−îc : 
⋅⋅⋅+++′=−= 6644 ha16
5
ha
4
1
)x(f
4
)1,1(D)1,2(D4
)2,2(D (11) 
Trong khi D(1,1) vµ D(2,1) sai kh¸c f′(x) phô thuéc vµo h2 th× D(2,2) sai kh¸c f′(x) phô 
thuéc vµo h4 . B©y giê ta l¹i chia ®«i b−íc h vµ nhËn ®−îc : 
D f x a h a h(2, ) ( ) ( / ) ( / ) ...2
1
4 2
5
16 24
4
6
6= + + +′ (12) 
vµ khö sè h¹ng cã h4 b»ng c¸ch t¹o ra : 
 D D f x a h(2, ) ( , ) ( ) ( ) ...2 16 32 15
15
64 6
6− − ′= + + + (13) 
Chia hai vÕ cña (13) cho -15 ta cã : 
D
D D
f x a h(3, )
(3, ) (2, )
( ) . ...3
16 2 2
15
1
64 6
6= = − −
− ′ (14) 
205
Víi lÇn tÝnh nµy sai sè cña ®¹o hµm chØ cßn phô thuéc vµo h6 . L¹i tiÕp tôc chia ®«i b−íc h 
vµ tÝnh D(4,4) th× sai sè phô thuéc h8 . S¬ ®å tÝnh ®¹o hµm theo ph−¬ng ph¸p Romberg lµ : 
 D(1,1) 
 D(2,1) D(2,2) 
 D(3,1) D(3,2) D(3,3) 
 D(4,1) D(4,2) D(4,3) D(4,4) 
 . . . . . . . . . . . . 
trong ®ã mçi gi¸ trÞ sau lµ gi¸ trÞ ngo¹i suy cña gi¸ trÞ tr−íc ®ã ë hµng trªn . 
Víi 2 ≤ j ≤ i ≤ n ta cã : 
D j
D j D jj
j(i, )
(i, ) (i , )
=
−
−
− − − −
−
1
1
4 1 1 1
4 1
vµ gi¸ trÞ khëi ®Çu lµ : 
D h h f x h f x hi i i i
(i, ) ( ) [ ( ) ( )]1 12= = + − −ϕ
víi hi = h/2
i-1 . 
Chóng ta ngõng l¹i khi hiÖu gi÷a hai lÇn ngo¹i suy ®¹t ®é chÝnh x¸c yªu cÇu. 
VÝ dô : T×m ®¹o hµm cña hµm f(x) = x2 + arctan(x) t¹i x = 2 víi b−íc tÝnh h = 0.5 . TrÞ chÝnh 
x¸c cña ®¹o hµm lµ 4.2 
201843569.4)]75.1(f)25.2(f[
25.02
1
)1,2(D
207496266.4)]5.1(f)5.2(f[
5.02
1
)1,1(D
=−×=
=−×= 
200458976.4)]875.1(f)125.2(f[
125.02
1
)1,3(D =−×= 
200492284.4
14
)2,2(D)2,3(D4
)3,3(D
200458976.4
14
)1,2(D)1,3(D4
)2,3(D
19995935.4
14
)1,1(D)1,2(D4
)2,2(D
21
2
1
1
1
1
==
==
==
−
−
−
−
−
−
Ch−¬ng tr×nh tÝnh ®¹o hµm nh− d−íi ®©y . Dïng ch−¬ng tr×nh tÝnh ®¹o hµm cña hµm 
cho trong function víi b−íc h = 0.25 t¹i xo = 0 ta nhËn ®−îc gi¸ trÞ ®¹o hµm lµ 1.000000001. 
Ch−¬ng tr×nh12-.1 
//Daoham_Romberg; 
#include 
#include 
#include 
#define max 11 
float h; 
void main() 
 { 
 float d[max]; 
 int j,k,n; 
 float x,p; 
 float y(float),dy(float); 
206
 clrscr(); 
 printf("Cho diem can tim dao ham x = "); 
 scanf("%f",&x); 
 printf("Tinh dao ham theo phuong phap Romberg\n"); 
 printf("cua ham f(x) = th(x) tai x = %4.2f\n",x); 
 n=10; 
 h=0.2; 
 d[0]=dy(x); 
 for (k=2;k<=n;k++) 
 { 
 h=h/2; 
 d[k]=dy(x); 
 p=1.0; 
 for (j=k-1;j>=1;j--) 
 { 
 p=4*p; 
 d[j]=(p*d[j+1]-d[j])/(p-1); 
 } 
 } 
 printf("y'= %10.5f\n",d[1]); 
 getch(); 
 } 
float y(float x) 
 { 
 float a=(exp(x)-exp(-x))/(exp(x)+exp(-x)); 
 return(a); 
 } 
float dy(float x) 
 { 
 float b=(y(x+h)-y(x-h))/(2*h); 
 return(b); 
 } 
§2. Kh¸i niÖm vÒ tÝch ph©n sè 
 Môc ®Ých cña tÝnh tÝch ph©n x¸c ®Þnh lµ ®¸nh gi¸ ®Þnh l−îng biÓu thøc : 
J f x
a
b
= ∫ ( )dx 
trong ®ã f(x) lµ hµm liªn tôc trong kho¶ng [a,b] vµ cã 
thÓ biÓu diÔn bëi ®−êng cong y= f(x). Nh− vËy tÝch 
ph©n x¸c ®Þnh J lµ diÖn tÝch SABba , giíi h¹n bëi ®−êng 
cong f(x) , trôc hoµnh , c¸c ®−êng th¼ng x = a vµ x = b 
. NÕu ta chia ®o¹n [a,b] thµnh n phÇn bëi c¸c ®iÓm xi 
th× J lµ gíi h¹n cña tæng diÖn tÝch c¸c h×nh ch÷ nhËt 
f(xi).(xi+1 - xi) khi sè ®iÓm chia tiÕn tíi ∝, nghÜa lµ : 
a a
b
A
B
y
x
207
J f x x x
n
i
i
n
i i
= −
→∞ = +
∑lim ( )( )
0
1
 NÕu c¸c ®iÓm chia xi c¸ch ®Òu , th× ( xi+1- xi ) = 
h . Khi ®Æt f(xo) = fo , f(x1) = f1 ,... ta cã tæng : 
n i
i
n
S h f=
=
∑
0
 Khi n rÊt lín , Sn tiÕn tíi J . Tuy nhiªn sai sè lµm trßn l¹i ®−îc tÝch luü . Do vËy cÇn 
ph¶i t×m ph−¬ng ph¸p tÝnh chÝnh x¸c h¬n . Do ®ã ng−êi ta Ýt khi dïng ph−¬ng ph¸p h×nh 
ch÷ nhËt nh− võa nªu . 
§3. Ph−¬ng ph¸p h×nh thang 
Trong ph−¬ng ph¸p h×nh thang , thay v× chia diÖn tÝch SABba thµnh c¸c h×nh ch÷ nhËt , 
ta l¹i dïng h×nh thang . VÝ dô nÕu chia thµnh 3 ®o¹n nh− h×nh vÏ th× : 
 S3 = t1 + t2 + t3 
trong ®ã ti lµ c¸c diÖn tÝch nguyªn tè . Mçi diÖn tÝch nµy lµ mét h×nh thang : 
 ti = [f(xi) + f(xi-1)]/ (2h) 
 = h(fi - fi-1) / 2 
Nh− vËy : 
 S3 = h[(fo+f1)+(f1+f2)+(f2+f3)] / 2 
 = h[fo+2f1+2f2+f3] / 2 
 Mét c¸ch tæng qu¸t chóng ta cã : 
)f2f2f2f(n
ab
S n1n1on ++⋅⋅⋅++= −
− 
hay : 
}f2ff{n
ab
S
n
1i
ion n ∑++− ==
 Mét c¸ch kh¸c ta cã thÓ viÕt : 
 f x dx f x hf a kh f a k h
a
b
a kh
a k h
k
n
k
n
( ) ( )dx { ( ) / [ ( ) ] / }
( )
∫ ∫∑ ∑= ≈ + + + +
+
+ +
=
−
=
−1
1
1
0
1
2 1 2 
hay : 
 f x h f a f a h f a n h f b
a
b
( )dx { ( ) / ( ) [ ( ) ] ( ) / }= + + + ⋅ ⋅ ⋅ + + − +∫ 2 1 2 
Ch−¬ng tr×nh tÝnh tÝch ph©n theo ph−¬ng ph¸p h×nh thang nh− sau : 
Ch−¬ng tr×nh 12-2 
//tinh tich phan bang phuong phap hinh_thang; 
#include 
#include 
#include 
float f(float x) 
 { 
 float a=exp(-x)*sin(x); 
 return(a); 
 }; 
208
void main() 
 { 
 int i,n; 
 float a,b,x,y,h,s,tp; 
 clrscr(); 
 printf("Tinh tich phan theo phuong phap hinh thang\n"); 
 printf("Cho can duoi a = "); 
 scanf("%f",&a); 
 printf("Cho can tren b = "); 
 scanf("%f",&b); 
 printf("Cho so buoc n = "); 
 scanf("%d",&n); 
 h=(b-a)/n; 
 x=a; 
 s=(f(a)+f(b))/2; 
 for (i=1;i<=n;i++) 
 { 
 x=x+h; 
 s=s+f(x); 
 } 
 tp=s*h; 
 printf("Gia tri cua tich phan la : %10.6f\n",tp); 
 getch(); 
 } 
 Dïng ch−¬ng tr×nh nµy tÝnh tÝch ph©n cña hµm cho trong function trong kho¶ng [0 , 
1] víi 20 ®iÓm chia ta cã J = 0.261084. 
§4. C«ng thøc Simpson 
 Kh¸c víi ph−¬ng ph¸p h×nh thang , ta chia ®o¹n [a,b] thµnh 2n phÇn ®Òu nhau bëi 
c¸c ®iÓm chia xi : 
 a = xo < x1 < x2 < ....< x2n = b 
 xi = a+ih ; h = (b - a)/ 2n víi i =0 , . . , 2n 
 Do yi = f(xi) nªn ta cã : 
∫∫∫ ∫
−
+++=
x
x
fdx...
x
x
fdx
b
a
x
x
fdxdx)x(f
n2
2n2
4
2
2
0
 §Ó tÝnh tÝch ph©n nµy ta thay hµm f(x) ë vÕ ph¶i b»ng ®a thøc néi suy Newton tiÕn 
bËc 2 : 
y
t2
)1t(t
yty)x(P 0
2
002 ∆−∆ ++= 
vµ víi tÝch ph©n thø nhÊt ta cã : 
dx)x(Pdx)x(f
x
x
x
x
2
0
2
0
2∫∫ = 
 §æi biÕn x = x0+th th× dx = hdt , víi x0 th× t =0 vµ víi x2 th× t = 2 nªn : 
209
|]y)
2
t
3
t(
2
1y
2
tty[h
dt)y
2
)1t(1
yty(hdx)x(P
2t
0t0
2
23
0
2
0
0
2
0
2
0
02
x
x
2
0
=
=∆∆
∆−∆∫∫
−++=
++= 
]yy4y[
3
h
]y)
2
4
3
8(
2
1y2y2[h
210
0
2
00
++=
−++= ∆∆ 
§èi víi c¸c tÝch ph©n sau ta còng cã kÕt qu¶ t−¬ng tù : 
]yy4y[
3
hdx)x(f 2i21i2i2
x
x
2i2
i2
++ ++=∫
+ 
Céng c¸c tÝch ph©n trªn ta cã : 
]y)yyy(2)yyy(4y[
3
hdx)x(f n22n2421n231o
b
a
++⋅⋅⋅++++⋅⋅⋅+++= −−∫ 
Ch−¬ng tr×nh dïng thuËt to¸n Simpson nh− sau : 
Ch−¬ng tr×nh 12-3 
//Phuong phap Simpson; 
#include 
#include 
#include 
float y(float x) 
 { 
 float a=4/(1+x*x); 
 return(a); 
 } 
void main() 
 { 
 int i,n; 
 float a,b,e,x,h,x2,y2,x4,y4,tp; 
 clrscr(); 
 printf("Tinh tich phan theo phuong phap Simpson\n"); 
 printf("Cho can duoi a = "); 
 scanf("%f",&a); 
 printf("Cho can tren b = "); 
 scanf("%f",&b); 
 printf("Cho so diem tinh n = "); 
 scanf("%d",&n); 
 h=(b-a)/n; 
 x2=a+h; 
 x4=a+h/2; 
 y4=y(x4); 
 y2=y(x2); 
 for (i=1;i<=n-2;i++) 
 { 
210
 x2+=h; 
 x4+=h; 
 y4+=y(x4); 
 y2+=y(x2); 
 } 
 y2=2*y2; 
 y4=4*(y4+y(x4+h)); 
 tp=h*(y4+y2+y(a)+y(b))/6; 
 printf("Gia tri cua tich phan la : %10.8f\n",tp); 
 getch(); 
 } 
 Dïng ch−¬ng tr×nh nµy tÝnh tÝch ph©n cña hµm trong function trong ®o¹n [0,1] víi 20 
kho¶ng chia cho ta kÕt qu¶ J = 3.14159265.