Giáo trình C++ - Chương 11: Nội suy và xấp xỉ hàm
Ch−¬ng 11 : néi suy vµ xÊp xØ hµm
§1.Néi suy Lagrange
Trong thùc tÕ nhiÒu khi ph¶i phôc håi mét hµm y = f(x) t¹i mäi gi¸ trÞ x trong mét
®o¹n [ a,b ] nµo ®ã mµ chØ biÕt mét sè nhÊt ®Þnh c¸c gi¸ trÞ cña hµm t¹i mét sè ®iÓm cho
tr−íc.C¸c gi¸ trÞ nµy ®−îc cung cÊp qua thùc nghiÖm hay tÝnh to¸n.V× vËy n¶y sinh vÊn ®Ò
to¸n häc lµ trªn ®o¹n a ≤ x ≤ b cho mét lo¹t c¸c ®iÓm xi ( i= 0,1,2..) vµ t¹i c¸c ®iÓm xi nµy
gi¸ trÞ cña hµm lµ yi = f(xi) ®· biÕt.B©y giê ta cÇn t×m ®a thøc :
Pn(x) = aoxn + a1xn-1 + …+an-1x + an
sao cho Pn(xi) = f(xi) = yi.§a thøc Pn(x) ®−îc gäi lµ ®a thøc néi suy cña hµm y = f(x).Ta chän
®a thøc ®Ó néi suy hµm y = f(x) v× ®a thøc lµ lo¹i hµm ®¬n gi¶n,lu«n cã ®¹o hµm vµ nguyªn
hµm.ViÖc tÝnh gi¸ trÞ cña nã theo thuËt to¸n Horner còng ®¬n gi¶n.
B©y giê ta x©y dùng ®a thøc néi suy kiÓu Lagrange.Gäi Li lµ ®a thøc :
(x − x0 )...(x − xi−1)(x − xi+1)...(x − xn )
(xi − x0 )...(xi − xi−1)(xi − xi+1 )...(xi − xn )
Li =
Râ rµng lµ Li(x) lµ mét ®a thøc bËc n vµ :
⎧⎪1
j=i
j≠i
=
⎨⎪0
(x )
Li
j
⎩
Ta gäi ®a thøc nµy lµ ®a thøc Lagrange c¬ b¶n.
B©y giê ta xÐt biÓu thøc :
n
=
(x)
f(x )L (x)
∑
Pn
i
i
i=0
Ta thÊy Pn(x) lµ mét ®a thøc bËc n v× c¸c Li(x) lµ c¸c ®a thøc bËc n vµ tho¶ m·n ®iÒu
kiÖn Pn(xi) = f(xi) = yi.Ta gäi nã lµ ®a thøc néi suy Lagrange.
Víi n = 1 ta cã b¶ng
x
y
x0
y0
x1
y1
§a thøc néi suy sÏ lµ :
P1(x) = yoL0(x) + y1L1(x1)
x − x1
x0 − x1
x − x0
L0 =
L1 =
x1 − x0
x − x0
x − x1
x0 − x1
nªn
P1(x) = y0
+ y
1 x1 − x0
Nh− vËy P1(x) lµ mét ®a thøc bËc nhÊt ®èi víi x
Víi n = 2 ta cã b¶ng
x
y
x0
y0
x1
y1
x2
y2
§a thøc néi suy sÏ lµ :
P2(x) = yoL0(x) + y1L1(x1) + y2L2(x2)
(x − x1 )(x − x2 )
(x0 − x1 )(x0 − x2 )
L0 =
(x − x0 )(x − x2 )
(x1 − x0 )(x1 − x2 )
L1 =
180
(x − x0 )(x − x1 )
(x2 − x0 )(x2 − x1)
L2 =
Nh− vËy P1(x) lµ mét ®a thøc bËc hai ®èi víi x
Trªn c¬ së thuËt to¸n trªn ta cã ch−¬ng tr×nh t×m ®a thøc néi suy cña mét hµm khi
cho tr−íc c¸c ®iÓm vµ sau ®ã tÝnh trÞ sè cña nã t¹i mét gi¸ trÞ nµo ®ã nh− sau :
Ch−¬ng tr×nh 11-1
#include <conio.h>
#include <stdio.h>
#include <ctype.h>
#define max 21
int maxkq,n;
float x[max],y[max],a[max],xx[max],yy[max];
float x0,p0;
void main()
{
int i,k;
char ok ;
void vaosolieu(void);
float lagrange(int,float [],float [],float);
void inkq(void);
clrscr();
printf("%24cNOI SUY DA THUC LAGRANGE\n",' ');
vaosolieu();
k=0;
ok='c';
while (ok=='c')
{
printf("Tinh gia tri cua y voi x la x0 = ");
scanf("%f",&x0);
p0=lagrange(n,x,y,x0);
printf("Gia tri cua y = %15.5f\n",p0);
printf("\n");
k=k+1;
maxkq=k;
xx[k]=x0;
yy[k]=p0;
flushall();
printf("Tinh tiep khong(c/k)?");
scanf("%c",&ok);
}
inkq();
181
}
void vaosolieu()
{
int i,t;
char ok;
printf("\n");
printf("Ham y = f(x)\n");
printf("So cap (x,y) nhieu nhat la max = 20\n");
printf("So diem da cho truoc n = ");
scanf("%d",&n);
for (i=1;i<=n;i++)
{
printf("x[%d] = ",i);
scanf("%f",&x[i]);
printf("y[%d] = ",i);
scanf("%f",&y[i]);
}
printf("\n");
printf(" SO LIEU BAN VUA NHAP\n");
printf("
x
y\n");
for (i=1;i<=n;i++)
printf("%8.4f
ok=' ';
%8.4f\n",x[i],y[i]);
t=1;
flushall();
while (t)
{
printf("\nCo sua so lieu khong(c/k):?");
scanf("%c",&ok);
if (toupper(ok)=='C')
{
printf("Chi so cua phan tu can sua i = ");
scanf("%d",&i);
printf("Gia tri moi : ");
printf("x[%d] = ",i);
scanf("%f",&x[i]);
printf("y[%d] = ",i);
scanf("%f",&y[i]);
flushall();
}
if (toupper(ok)!='C')
t=0;
}
}
float lagrange(int n,float x[max],float y[max],float x0)
{
int i,k;
182
float g0;
p0=0.0;
for (k=1;k<=n;k++)
{
g0=1.0;
for (i=1;i<=n;i++)
if (i!=k)
g0=g0*(x0-x[i])/(x[k]-x[i]);
p0=p0+y[k]*g0;
}
return(p0);
}
void inkq()
{
int i,j,k;
printf("\n");
printf("%24cBANG SO LIEU\n",' ');
printf("%18cx %24cy\n",' ',' ');
for (i=1;i<=n;i++)
printf("%20.4f %25.4f\n",x[i],y[i]);
printf("\n");
printf("%24cKET QUA TINH TOAN\n",' ');
printf("%14cx %10cy\n",' ',' ');
for (k=1;k<=maxkq;k++)
printf("%15.5f %15.5f\n",xx[k],yy[k]);
getch();
}
Gi¶ sö ta cã b¶ng c¸c gi¸ trÞ x,y :
x
y
0
0
3
-2
2
-2
4
4
-3.75 10
vËy theo ch−¬ng tr×nh t¹i x = 2.5 y = -3.3549.
§2.Néi suy Newton
B©y giê ta xÐt mét c¸ch kh¸c ®Ó x©y dùng ®a thøc néi suy gäi lµ ph−¬ng ph¸p
Newton.Tr−íc hÕt ts ®−a vµo mét kh¸i niÖm míi lµ tØ hiÖu
Gi¶ sö hµm y = y(x) cã gi¸ trÞ cho trong b¶ng sau :
x
y
x0
y0
x1
y1
x2
y2
xn-1
yn-1
xn
yn
…
…
TØ hiÖu cÊp 1 cña y t¹i xi,xj lµ :
183
yi − yj
xi − xj
y[xi ,xj ] =
TØ hiÖu cÊp hai cña y t¹i xi,xj,xk lµ :
y[xi ,xj ]− y[xj,xk ]
y[xi ,xj,xk ] =
xi − xk
v.v.
Víi y(x) = Pn(x) lµ mét ®a thøc bËc n th× tØ hiÖu cÊp 1 t¹i x,x0 :
Pn (x) − Pn (x0 )
Pn [x,x0 ] =
x − x0
lµ mét ®a thøc bËc (n-1).TØ hiÖu cÊp 2 t¹i x,x0,x1 :
Pn[x,x0 ]− Pn [x0 ,x1]
Pn [x,x0 ,x1] =
x − x1
lµ mét ®a thøc bËc (n-2) v.v vµ tíi tØ hiÖu cÊp (n+1) th× :
Pn[ x,xo,..,xn] = 0
Tõ c¸c ®Þnh nghÜa tØ hiÖu ta suy ra :
Pn(x) = Pn(x0) + ( x- x0)Pn[x,xo]
Pn[x,x0] = Pn[x0,x1] + ( x- x1) Pn[x,xo,x1]
Pn[x,xo,x1] = Pn[x0,x1,x2] + ( x- x2) Pn[x,xo,x1,x2]
............
Pn[x,xo,..,xn-1] = Pn[x0,x1,..,xn] + ( x- xn) Pn[x,xo,..,xn]
Do Pn[ x,xo,..,xn] = 0 nªn tõ ®ã ta cã :
Pn(x) = Pn(x0) + (x - x0)Pn[xo,x1] + (x - x0)(x - x1)Pn[x0,x1,x2] +…
+(x - x0)…(x - xn-1)Pn[x0,…,xn]
NÕu Pn(x) lµ ®a thøc néi suy cña hµm y=f(x) th× :
Pn(xi) = f(xi) = yi víi i = 0 ÷ n
Do ®ã c¸c tØ hiÖu tõ cÊp 1 ®Õn cÊp n cña Pn vµ cña y lµ trïng nhau vµ nh− vËy ta cã :
Pn(x) = y0 + (x - x0)y[x0,x1] + (x - x0)(x - x1)y[x0,x1,x2] +..+
(x - x0)(x - x1)...(x - xn-1)y[x0,..,xn]
§a thøc nµy gäi lµ ®a thøc néi suy Newton tiÕn xuÊt ph¸t tõ nót x0 cña hµm y =
f(x).Ngoµi ®a thøc tiÕn cßn cã ®a thøc néi suy Newton lïi xuÊt ph¸t tõ ®iÓm xn cã d¹ng nh−
sau :
Pn(x) = yn + (x - xn)y[xn,xn-1] + (x - xn)(x - xn-1)y[xn,xn-1,xn-2] +..+
(x - xn)(x - xn-1)...(x - x1)y[xn,..,x0]
Tr−êng hîp c¸c nót c¸ch ®Òu th× xi = x0 +ih víi i = 0,1,..,n.Ta gäi sai ph©n tiÕn cÊp 1
t¹i i lµ :
∆yi = yi+1 - yi
vµ sai ph©n tiÕn cÊp hai t¹i i :
∆2yi = ∆(∆yi) = yi+2 - 2yi+1 + yi
.........
vµ sai ph©n tiÕn cÊp n lµ :
∆nyi = ∆(∆n-1yi)
Khi ®ã ta cã :
= ∆y0
y[x
0,x1]
h
2 y
∆
0
=
y[x
0,x1,x2 ]
2h2
...........
184
n y
∆
0
= (n!hn )
y[x ,...,x
]
0
n
B©y giê ®Æt x = x0 + ht trong ®a thøc Newton tiÕn ta ®−îc :
t(t−1)...(t−n+1)
+ t(t−1)
2 y
n y
∆
+
=
+
+ ... +
y0 t∆y0
(x ht)
Pn
∆
0
0
0
2!
n!
th× ta nhËn ®−îc ®a thøc Newton tiÕn xuÊt ph¸t tõ x0 trong tr−êng hîp nót c¸ch ®Òu.Víi n =1
ta cã :
P1(x0+ht) = y0 + ∆y0
Víi n =2 ta cã :
+ t(t−1)
+
=
+
y0 t∆y0
2 y
∆
(x ht)
P2
0
0
2
Mét c¸ch t−¬ng tù ta cã kh¸i niÖm c¸c sai ph©n lïi t¹i i :
∇yi = yi - yi-1
∇2yi = ∇(∇yi) = yi - 2yi-1 + yi-2
.........
∇nyi = ∇(∇n-1yi)
vµ ®a thøc néi suy Newton lïi khi c¸c ®iÓm néi suy c¸ch ®Òu :
t(t+1)...(t+n−1)
+ t(t+1)
2 y
n y
∇
+
=
+
+ ... +
yn t∇yn
(x ht)
Pn
∇
0
n
n
2!
n!
VÝ dô : Cho hµm nh− b¶ng sau :
x
y
0.1
0.2
0.3
0.4
0.09983 0.19867 0.29552 0.38942
Ta tÝnh gi¸ trÞ cña hµm t¹i 0.14 b»ng ®a thøc néi suy Newton v× c¸c mèc c¸ch ®Òu
h = 0.1.Ta cã b¶ng sai ph©n sau :
∆2y
∆3y
i
0
x
0.1
y
∆y
0.09983
0.09884
1
2
3
0.2
0.3
0.4
0.19867
0.29552
0.38942
-
0.00199
0.09685
0.09390
-0.00096
-
0.00295
Ta dïng c«ng thøc Newton tiÕn víi ®iÓm gèc lµ x0 = 0.1.h = 0.1.Víi x = 0.14 ta cã
0.14 = 0.1 + 0.1t nªn t = 0.4 vµ kÕt qu¶ lµ :
t(t −1)
t(t −1)(t − 2)
P(0.1+ 0.1t) =0.09983+ t.0.099884+
tr×nh néi suy Newton nh− sau :
0.00199−
0.00096=0.13954336Ch−¬ng
2!
3!
Ch−¬ng tr×nh 11-2
185
//Noi suy Newton
#include <conio.h>
#include <stdio.h>
#include <ctype.h>
#define max 11
void main()
{
int i,j,k,n,t;
float a[max],b[max],x[max],y[max];
char ok;
float x0,p;
clrscr();
printf("So diem da cho n = ");
scanf("%d",&n);
for (i=1;i<=n;i++)
{
printf("x[%d] = ",i);
scanf("%f",&x[i]);
printf("y[%d] = ",i);
scanf("%f",&y[i]);
}
printf("%10cBANG SO LIEU\n",' ');
printf("%8cx%30cy\n",' ',' ');
for (i=1;i<=n;i++)
printf("%4c%8.4f%23c%8.4f\n",' ',x[i],' ',y[i]);
ok=' ';
t=0;
flushall();
while (t)
{
printf("Co sua so lieu khong(c/k): ");
scanf("%c",&ok);
if (toupper(ok)=='C')
{
printf("Chi so cua phan tu can sua i = ");
scanf("%d",&i);
printf("Gia tri moi : ");
printf("x[%d] = ",i);
scanf("%f",&x[i]);
printf("y[%d] = ",i);
scanf("%f",&y[i]);
flushall();
}
if (toupper(ok)!='C')
t=0;
}
a[1]=y[1];
for (j=1;j<=n-1;j++)
186
{
}
for (i=1;i<=n-j;i++)
y[i]=(y[i+1]-y[i])/(x[i+j]-x[i]);
a[j+1]=y[1];
b[n]=a[n];
for (k=n-1;k>=1;k--)
{
for (j=n-1;j>=1;j--)
b[j]=a[j] ;
for (i=n-1;i>=k;i--)
a[i]=a[i]-b[i+1]*x[k];
}
for (i=n;i>=1;i--)
printf("He so bac %d la :%8.4f\n",i-1,a[i]);
printf("\n");
k=0;
ok='c';
flushall();
while (ok=='c')
{
printf("Tinh gia tri cua y tai x = ");
scanf("%f",&x0);
p=0;
for (k=n;k>=1;k--)
p=p*x0+a[k];
printf("Tri so noi suy tai x0 = %4.2f la : %10.5f\n",x0,p);
getch();
printf("Ban co muon tinh tiep cac diem khac khong(c/k)");
do
scanf("%c",&ok);
while ((ok!='c')&&(ok!='k'));
}
}
Dïng ch−¬ng tr×nh nµy néi suy c¸c gi¸ trÞ cho trong b¶ng sau
0
1
0.2
0.4
0.6
0.8
1.0
1.2214027 1.4918247 1.8221188 2.2255409 2.7182818
6
3
3
ta cã c¸c hÖ sè cña ®a thøc néi suy : 0.0139(bËc 5),0.0349(bËc 4),0.1704(bËc3),0.4991(bËc
2),1.0001(bËc 1) vµ 1.0000(bËc 0).
§3.Néi suy Aitken
Mét d¹ng kh¸c cña ®a thøc néi suy ®−îc x¸c ®Þnh b»ng thuËt to¸n Aitken.Gi¶ sö ta cã
n ®iÓm ®· cho cña hµm f(x).Nh− vËy qua hai ®iÓm x0 vµ x1 ta cã ®a thøc néi suy Lagrange
cña hµm f(x) ®−îc viÕt d−íi d¹ng :
187
y0 x0 − x
y1 x1 − x
P01(x) =
x1 − x0
lµ mét ®a thøc bËc 1 :
x − x1
x0 − x1
x − x0
P01(x) = y0
+ y
1 x1 − x0
.Khi x = x0 th× :
y0 x0 − x0
y1 x1 − x0
P01(x0 ) =
= y0
x1 − x0
Khi x = x1 th× :
y0 x0 − x1
y1 x1 − x1
P01(x1 ) =
= y1
x1 − x0
§a thøc néi suy Lagrange cña f(x) qua 3 ®iÓm x0,x1,x2 cã d¹ng :
P01(x) x0 − x
P12 (x) x2 − x
P012 (x) =
x2 − x0
vµ lµ mét ®a thøc bËc 2:
(x − x1)(x − x2 )
(x0 − x1)(x0 − x2 )
(x − x0 )(x − x2 )
(x − x0 )(x − x1 )
(x2 − x0 )(x2 − x1)
P012 (x) = y0
+ y
+ y2
1 (x1 − x0 )(x1 − x2 )
Khi x = x0 th× :
y0
x0 − x0
P12 (x) x2 − x0
P012 (x0 ) =
= y0
x2 − x0
Khi x = x1 th× :
y1 x0 − x1
y1 x2 − x1
P012 (x1 ) =
= y1
x2 − x0
Khi x = x2 th× :
P01(x2 ) x0 − x2
y2 x2 − x2
x2 − x0
P012 (x2 ) =
= y2
Tæng qu¸t ®a thøc néi suy Lagrange qua n ®iÓm lµ :
P01..(n−1) (x) x0 − x
P12..n (x)
xn − x
P012..n (x) =
x2 − x0
Nh− vËy ta cã thÓ dïng phÐp lÆp ®Ó x¸c ®Þnh lÇn l−ît c¸c ®a thøc Lagrange.S¬ ®å tÝnh
to¸n nh− vËy gäi lµ s¬ ®å Neville-Aitken.
VÝ dô : Cho c¸c cÆp ®iÓm (0,0.4),(1.4,1.5),(2.6,1.8),(3.9,2.6),tÝnh y t¹i x=2
188
y0 x0 − x
y1 x1 − x
0.4 − 2
1.5 − 0.6
P01(x) =
P12 (x) =
P012 (x) =
P23 (x) =
=
=
= 1.97143
= 1.65
x1 − x0
1.4 − 0
y1 x1 − x
y2 x2 − x
1.5 − 0.6
1.8 0.6
x2 − x1
2.6 −1.4
P01(x) x0 − x
P12 (x) x2 − x
1.97143 − 2
1.65 0.6
2.6 − 0
=
= 1.7242
x2 − x0
y2 x2 − x
y3 x3 − x
1.8 0.6
2.6 1.9
=
= 1.4308
x3 − x2
3.9 − 2.6
P12 (x) x1 − x
P23 (x) x3 − x
1.65 − 0.6
1.4308 1.9
P123 (x) =
=
= 1.5974
x3 − x1
3.9 −1.4
P012 (x) x0 − x
P123 (x) x3 − x
1.7242 − 2
1.5974 1.9
P0123 (x) =
=
= 1.6592
x3 − x0
3.9 − 0
Ch−¬ng tr×nh ®−îc viÕt nh− sau
Ch−¬ng tr×nh 11-3
//Noi suy Aitken
#include <conio.h>
#include <stdio.h>
#include <ctype.h>
#define max 11
void main()
{
float x[max],y[max],yd[max];
float x1;
int j,k,n,n1;
clrscr();
printf("Cho so diem da co n = ");
scanf("%d",&n1);
n=n1-1 ;
for (k=0;k<=n;k++)
{
printf("x[%d] = ",k+1);
scanf("%f",&x[k]);
printf("y[%d] = ",k+1);
scanf("%f",&y[k]);
}
printf("Cho diem can tinh gia tri cua ham x1 = ");
189
scanf("%f",&x1);
for (k=0;k<=n-1;k++)
{
yd[k]=(y[k]*(x1-x[k+1])-y[k+1]*(x1-x[k]))/(x[k]-x[k+1]);
if (k!=0)
for (j=k-1;j>=0;j--)
yd[j]=(yd[j]*(x1-x[k+1])-yd[j+1]*(x1-x[j]))/(x[j]-x[k+1]);
}
printf("Gia tri ham tai x = %6.3f la y = %8.4f\n",x1,yd[0]);
getch();
}
Dïng ch−¬ng tr×nh nµy ®Ó néi suy c¸c cÆp sè (1,3),(2,5),(3,7),(4,9) vµ (5,11) t¹i x =
2.5 ta cã y = 6.
§4.XÊp xØ hµm b»ng ph−¬ng ph¸p b×nh ph−¬ng bÐ nhÊt
Trong c¸c môc tr−íc ta ®· néi suy gi¸ trÞ cña hµm.Bµi to¸n ®ã lµ cho mét hµm d−íi
d¹ng b¶ng sè vµ ph¶i t×m gi¸ trÞ cña hµm t¹i mét gi¸ trÞ cña ®èi sè kh«ng n»m trong b¶ng.
Trong thùc tÕ,bªn c¹nh bµi to¸n néi suy ta cßn gÆp mét d¹ng bµi to¸n kh¸c.§ã lµ t×m
c«ng thøc thùc nghiÖm cña mét hµm.Néi dung bµi to¸n lµ tõ mét lo¹t c¸c ®iÓm cho tr−íc (cã
thÓ lµ c¸c gi¸ trÞ cña mét phÐp ®o nµo ®ã) ta ph¶i t×m mét hµm xÊp xØ c¸c gi¸ trÞ ®· cho.Ta sÏ
dïng ph−¬ng ph¸p b×nh ph−¬ng tèi thiÓu ®Ó gi¶i bµi to¸n.Gi¶ sö cã mÉu quan s¸t (xi,yi ) cña
hµm y = f(x).Ta chän hµm f(x) cã d¹ng :
f(x) = a0f0(x) + a1f1(x) + a2f2(x)...
(1)
Trong ®ã c¸c hµm f0(x),f1(x),f2(x) v.v.lµ (m+1) hµm ®éc lËp tuyÕn tÝnh mµ ta cã thÓ chän tuú
ý vµ c¸c hÖ sè ai lµ tham sè ch−a biÕt mµ ta ph¶i x¸c ®Þnh dùa vµo hÖ hµm ®· chän vµ c¸c
®iÓm quan s¸t.Sai sè gi÷a trÞ ®o ®−îc vµ trÞ tÝnh theo (1) lµ :
ei = yi - f(xi)
(2)
Sai sè nµy cã thÓ ©m hay d−¬ng tuú tõng gi¸ trÞ cña yi.Khi dïng ph−¬ng ph¸p b×nh ph−¬ng
bÐ nhÊt ta xÐt b×nh ph−¬ng cña sai sè t¹i mét ®iÓm :
e2i =
[
yi − f(xi ) 2
]
(3)
Víi n ®iÓm tæng b×nh ph−¬ng cña sai sè sÏ lµ :
n
n
2
S = e2 =
{
y −
[
a0f0 (xi ) + a1f1 (xi ) + ⋅⋅⋅ + anfn (xi )
]
}
∑ ∑
i
i
i=1
i=1
Râ rµng S lµ hµm cña c¸c gi¸ trÞ cÇn t×m ai.vµ chóng ta sÏ chän c¸c ai sao cho S ®¹t gi¸ trÞ
∂S
∂ai
min,nghÜa lµ c¸c ®¹o hµm
ph¶i b»ng kh«ng.Ta sÏ xÐt c¸c tr−êng hîp cô thÓ.
1.Hµm xÊp xØ cã d¹ng ®a thøc : Trong tr−êng hîp tæng qu¸t ta chän hÖ hµm xÊp xØ lµ mét
®a thøc,nghÜa lµ :
f(x) = a0 + a1x + a2x2 +...+ amxm
VËy hµm S lµ :
2
S =
(
yi −a0 + a1x + a2x +⋅⋅⋅+ amx
)
∂S
∂ai
Theo ®iÒu kiÖn ®¹o hµm
= 0ta nhËn ®−îc hÖ ph−¬ng tr×nh:
190
n
n
n
⎧
a
xm + a
x
m−1 + ⋅⋅⋅ +na =
y
i
⎪ ∑
∑
∑
m
i
m−1
i
0
i=1
i=1
i=1
⎪
n
n
n
n
m+1 + am−1 xm + ⋅⋅⋅ +a
x = x y
i i
i
⎪
am
∑
⎪
x
∑
0 ∑ ∑
i=1 i=1
n n
i
i
i=1
i=1
⎪
n
n
am xm+2 + am−1
x
m+1 + ⋅⋅⋅ +a0 x2 = x2y
⎪
∑
∑
∑ ∑
i
i
i
i
i
i
⎨
⎪
⎪
⎪
i=1
i=1
i=1
i=1
n
n
n
n
am xm+3 + am−1 xm+2 + ⋅⋅⋅ +a0 x3 = x3y
∑
∑
∑ ∑
i=1 i=1
i
i
i
i
i=1
i=1
⋅⋅⋅
⎪
n
n
n
n
⎪
a
⎪
⎩
x2m + am−1
x
2m−1 + ⋅⋅⋅ +a0 xm = xm y
∑
∑
∑
∑
i=1
m
i
i
i
i
i
i=1
i=1
i=1
§©y lµ mét hÖ ph−¬ng tr×nh tuyÕn tÝnh.Gi¶i nã ta nhËn ®−îc c¸c gÝa trÞ ai.Sau ®©y lµ
ch−¬ng tr×nh viÕt theo thuËt to¸n trªn.
Ch−¬ng tr×nh 11-4
//Xap xi da thuc
#include <conio.h>
#include <stdio.h>
#include <ctype.h>
#define max 11
void main()
{
int i,j,k,m,n,p,kp,t;
float a[max],x[max],y[max],y1[max];
float b[max][max];
char ok;
float s,sx,s1,c,d;
clrscr();
printf("PHUONG PHAP BINH PHUONG TOI THIEU");
printf("\n");
printf("Cho bac cua da thuc xap xi m = ");
scanf("%d",&m);
printf("So diem da cho n = ");
scanf("%d",&n);
for (i=1;i<=n;i++)
{
printf("x[%d] = ",i);
scanf("%f",&x[i]);
printf("y[%d] = ",i);
scanf("%f",&y[i]);
}
x[0]=1;
printf("\n");
printf("%4cBANG SO LIEU\n",' ');
191
printf("%8cx%30cy\n",' ',' ');
for (i=1;i<=n;i++)
printf("%4c%8.4f%20c%8.4f\n",' ',x[i],' ',y[i]);
ok=' ';
t=1;
flushall();
while (t)
{
printf("Co sua so lieu khong(c/k): ");
scanf("%c",&ok);
if (toupper(ok)=='C')
{
printf("Chi so cua phan tu can sua i = ");
scanf("%d",&i);
printf("Gia tri moi : ");
printf("x[%d] = ",i);
scanf("%f",&x[i]);
printf("y[%d] = ",i);
scanf("%f",&y[i]);
flushall();
}
if (toupper(ok)!='C')
t=0;
}
//for (i=0;i<=n;i++)
//a[i]=0.0;
printf("\n");
printf("CAC GIA TRI DA CHO");
printf("\n");
printf("X = ");
for (i=1;i<=n;i++)
printf("%c%8.3f",' ',x[i]);
printf("\n");
printf("Y = ");
for (i=1;i<=n;i++)
printf("%c%8.3f",' ',y[i]);
printf("\n");
for (p=0;p<=m;p++)
{
y1[p]=0.0;
for (i=1;i<=n;i++)
{
sx=1.0;
for (j=1;j<=p;j++)
sx*=x[i];
y1[p]+=y[i]*sx;
}
}
for (p=0;p<=m;p++)
for (k=0;k<=m;k++)
192
{
kp=k+p;
b[p][k]=0.0;
for (i=1;i<=n;i++)
{
sx=1.0;
for (j=1;j<=kp;j++)
sx*=x[i];
b[p][k]+=sx;
}
}
for (i=0;i<=m-1;i++)
{
c=1.0/b[i][i];
for (k=i+1;k<=m;k++)
{
d=b[i][k];
for (j=i+1;j<=m;j++)
b[k][j]-=b[i][j]*c*d;
y1[k]-=y1[i]*c*d;
b[i][k]*=c;
}
y1[i]*=c;
}
y1[m]/=b[m][m];
for (i=m-1;i>=0;i--)
for (j=i+1;j<=m;j++)
y1[i]-=b[i][j]*y1[j];
printf("\n");
printf("CAC HE SO CUA DA THUC CAN TIM");
printf("\n");
for (i=0;i<=m;i++)
printf("a[%d] = %10.5f\n",i,y1[i]);
getch();
}
Víi c¸c gi¸ trÞ x,y ®o ®−îc theo b¶ng
x
y
7
8
9
9,1
10
9,4
11
9,5
12
9,5
13
9,4
7,4 8,4
ta cã n = 7 vµ chän m = 2 vµ tÝnh ®−îc theo ch−¬ng tr×nh c¸c hÖ sè :
a0 = -0.111905 ; a1 = 2.545238 ; a2 = -4.857143
vµ hµm xÊp xØ sÏ lµ : f(x) = -0.111905 + 2.545238x -4.857143x2
2.Hµm d¹ng Aecx : Khi c¸c sè liÖu thÓ hiÖn mét sù biÕn ®æi ®¬n ®iÖu ta dïng hµm xÊp xØ lµ
y = Aecx.LÊy logarit hai vÕ ta cã :
lny = lnA + cxlne
193
∂S
∂ai
Theo ®iÒu kiÖn ®¹o hµm
= 0ta cã hÖ ph−¬ng tr×nh :
n
n
⎧
c
xi + n ln A = ln yi
⎪ ∑
∑
i=1
n
⎪
i=1
⎨
n
n
⎪
c
x2 + ln A x = x ln y
∑
∑ ∑
i
i
i
i
⎪
⎩ i=1
i=1
i=1
Gi¶i hÖ ph−¬ng tr×nh nµy ta cã c¸c hÖ sè A vµ c :
Ch−¬ng tr×nh 11-5
//xap_xi_e_mu;
#include <conio.h>
#include <stdio.h>
#include <ctype.h>
#include <math.h>
#define max 11
void main()
{
int i,n,t;
float x[max],y[max];
char ok;
float a,b,c,d,e,f,d1,d2,d3;
clrscr();
printf("PHUONG PHAP BINH PHUONG TOI THIEU");
printf("\n");
printf("So diem da cho n = ");
scanf("%d",&n);
for (i=1;i<=n;i++)
{
printf("x[%d] = ",i);
scanf("%f",&x[i]);
printf("y[%d] = ",i);
scanf("%f",&y[i]);
}
x[0]=1.0;
printf("%4cBANG SO LIEU\n",' ');
printf("%8cx%30cy\n",' ',' ');
for (i=1;i<=n;i++)
printf("%4c%8.4f%23c%8.4f\n",' ',x[i],' ',y[i]);
ok=' ';
t=1;
while (t)
{
printf("Co sua so lieu khong(c/k): ");
scanf("%c",&ok);
if (toupper(ok)=='C')
{
194
printf("Chi so cua phan tu can sua i = ");
scanf("%d",&i);
printf("Gia tri moi : ");
printf("x[%d] = ",i);
scanf("%f",&x[i]);
printf("y[%d] = ",i);
scanf("%f",&y[i]);
}
if (toupper(ok)!='C')
t=0;
}
printf("CAC GIA TRI DA CHO");
printf("\n");
printf("X = ");
for (i=1;i<=n;i++)
printf("%c%8.3f",' ',x[i]);
printf("\n");
printf("Y = ");
for (i=1;i<=n;i++)
printf("%c%8.3f",' ',y[i]);
printf("\n");
a=0.0;
for (i=1;i<=n;i++)
a+=x[i];
b=n;
c=0.0;
for (i=1;i<=n;i++)
c+=log(y[i]);
d=0.0;
for (i=1;i<=n;i++)
d+=x[i]*x[i];
e=0.0;
for (i=1;i<=n;i++)
e+=x[i]*log(y[i]);
d1=a*a-d*b;
d2=c*a-e*b;
d3=a*e-c*d;
c=d2/d1;
a=d3/d1;
printf("\n");
printf("He so A = %8.4f",exp(a));
printf(" va so mu c = %8.4",c);
printf("\n");
printf("\nBANG CAC GIA TRI TINH TOAN");
printf("\n");
printf("%5cx%28cy\n",' ',' ');
for (i=1;i<=n;i++)
{
printf("%8.4f%21c%8.4f\n",x[i],' ',exp(a)*exp(c*x[i]));
}
195
getch();
}
Víi c¸c gi¸ trÞ x,y ®o ®−îc theo b¶ng
x
0
2
4
6
162
8
76
10
43
12
19
y 128 635 324
0
ta cã n = 7 vµ tÝnh ®−îc theo ch−¬ng tr×nh c¸c hÖ sè : A = 1285.44 va c = -0.3476 vµ hµm
xÊp xØ sÏ lµ : f(x) = 1285.44
3.Hµm d¹ng Axq : Khi c¸c sè liÖu thÓ hiÖn mét sù biÕn ®æi ®¬n ®iÖu ta còng cã thÓ dïng
hµm xÊp xØ lµ y = Axq.LÊy logarit hai vÕ ta cã :
lny = lnA + qlnx
Theo ®iÒu kiÖn ®¹o hµm triÖt tiªu ta cã hÖ ph−¬ng tr×nh :
n
n
⎧
q
lnxi + n ln A = ln yi
⎪ ∑
∑
⎪
i=1
i=1
⎨
n
n
n
⎪
q
ln2 x + ln A lnx = lnx ln y
∑
∑
∑
i
i
i
i
⎪
⎩
i=1
i=1
i=1
Gi¶i hÖ ph−¬ng tr×nh nµy ta cã c¸c hÖ sè A vµ q :
Ch−¬ng tr×nh 11-6
//xap_xi_x_mu;
#include <conio.h>
#include <stdio.h>
#include <ctype.h>
#include <math.h>
#define max 11
void main()
{
int i,n,t;
float x[max],y[max];
char ok;
float a,b,c,d,e,f,d1,d2,d3;
clrscr();
printf("PHUONG PHAP BINH PHUONG TOI THIEU");
printf("\n");
printf("So diem da cho n = ");
scanf("%d",&n);
for (i=1;i<=n;i++)
{
printf("x[%d] = ",i);
scanf("%f",&x[i]);
printf("y[%d] = ",i);
scanf("%f",&y[i]);
196
}
x[0]=1.0;
printf("%4cBANG SO LIEU\n",' ');
printf("%8cx%30cy\n",' ',' ');
for (i=1;i<=n;i++)
printf("%4c%8.4f%23c%8.4f\n",' ',x[i],' ',y[i]);
ok=' ';
flushall();
t=1;
while (t)
{
printf("Co sua so lieu khong(c/k): ");
scanf("%c",&ok);
if (toupper(ok)=='C')
{
printf("Chi so cua phan tu can sua i = ");
scanf("%d",&i);
printf("Gia tri moi : ");
printf("x[",i,"] = ");
scanf("%f",&x[i]);
printf("y[%d] = ",i);
scanf("%f",&y[i]);
}
if (toupper(ok)!='C')
t=0;
}
printf("\n");
printf("\nCAC GIA TRI DA CHO");
printf("\n");
printf("X = ");
for (i=1;i<=n;i++)
printf("%c%8.3f",' ',x[i]);
printf("\n");
printf("Y = ");
for (i=1;i<=n;i++)
printf("%c%8.3f",' ',y[i]);
printf("\n");
a=0.0;
for (i=1;i<=n;i++)
a+=log(x[i]);
b=n;
c=0.0;
for (i=1;i<=n;i++)
c+=log(y[i]);
d=0.0;
for (i=1;i<=n;i++)
d+=log(x[i])*log(x[i]);
e=0.;
for (i=1;i<=n;i++)
e+=log(x[i])*log(y[i]);
197
d1=a*a-d*b;
d2=c*a-e*b;
d3=a*e-c*d;
c=d2/d1;
a=d3/d1;
printf("\n");
printf("He so A = %8.4f",exp(a));
printf(" va so mu q = %8.4f\n",c);
printf("\n");
printf("\nBANG CAC GIA TRI TINH TOAN\n");
printf("%5cx%27cy\n",' ',' ');
for (i=1;i<=n;i++)
{
printf("%8.4f%20c%8.4f\n",x[i],' ',exp(a)*exp(c*log(x[i])));
}
getch();
}
Víi c¸c gi¸ trÞ x,y ®o ®−îc theo b¶ng
x
1
2
4
5
6
y
7.1 27.8 62.1
110
161
ta cã n = 5 vµ tÝnh ®−îc theo ch−¬ng tr×nh c¸c hÖ sè : A = 7.1641 vµ q = 1.9531 vµ hµm xÊp
xØ sÏ lµ : f(x) = 1285.44x1.9531
4.Hµm l−îng gi¸c : Khi quan hÖ y=f(x) cã d¹ng tuÇn hoµn ta dïng hµm xÊp xØ lµ tæ hîp
tuyÕn tÝnh cña c¸c hµm sin vµ cosin d¹ng :
n
n
f(x) = a + a cos(iωx) + b sin(iωx)
∑
∑
0
i
i
i=1
i=1
§Ó ®¬n gi¶n tr−íc hÕt ta xÐt hµm chØ cã mét sè h¹ng sin-cos,nghÜa lµ :
f(x) = a0 + a1 cosωx + b1 sinωx
Hµm S sÏ cã d¹ng :
n
S =
[
y − (a + a cosωx + b sinωx) 2
]
∑
i
0
1
1
i=1
Theo ®iÒu kiÖn ®¹o hµm triÖt tiªu ta cã hÖ ph−¬ng tr×nh ®èi víi c¸c hÖ sè d¹ng :
⎡
⎢
⎢
⎤
⎥
n
cosωx
sinωx
⎡
⎢
⎢
⎤
⎥
⎥
y
∑
∑
a
∑
⎡ ⎤
0
cosωx
sinωx
cos2 ωx
cosωxsinωx a
=
ycosωx
ysinωx
⎢ ⎥
∑
∑
∑
∑
∑
∑
1
⎥
⎢ ⎥
cosωxsinωx
sin2 ωx
b1
⎢
⎣
⎥
⎣ ⎦
⎢
⎣
⎥
⎦
∑
∑
⎦
Do :
sinωx
cosωx
∑
∑
∑
∑
= 0
= 0
n
n
sin2 ωx
cos2 ωx
1
2
1
2
∑
=
=
n
n
cosωxsinωx
= 0
n
nªn hÖ ph−¬ng tr×nh cã d¹ng ®¬n gi¶n :
198
⎡
⎢
⎢
⎤
⎥
⎥
y
a
∑
⎡ ⎤
⎡
⎢
⎤
⎥
n
0
0
0
n 2
0
⎢ ⎥
0 n 2
a1
=
ycosωx
ysinωx
∑
∑
⎢ ⎥
b1
⎢
⎣
⎥
⎦
0
0
⎢
⎣
⎥
⎦
⎣ ⎦
Gi¶i hÖ ta cã :
y
2
2
∑
a0 =
a1 =
ycosωx
b1 =
ysinωx
∑
∑
n
n
n
Trong tr−êng hîp tæng qu¸t,mét c¸ch t−¬ng tù ta cã :
y
2
2
∑
a0 =
ai =
ycosiωx
bi =
ysiniωx
∑
∑
n
n
n
Ch−¬ng tr×nh t×m c¸c hÖ sè ai vµ bi ®−îc thÓ hiÖn nh− sau :
Ch−¬ng tr×nh 11-7
//xap_xi_sin_cos;
#include <conio.h>
#include <stdio.h>
#include <ctype.h>
#include <math.h>
#define max 11
#define pi 3.15159
void main()
{
int i,j,m,n,t;
float x[max],y[max],a[max],b[max];
char ok;
float omg,t1;
clrscr();
printf("PHUONG PHAP BINH PHUONG TOI THIEU");
printf("\n");
printf("Cho so so hang sin-cos m = ");
scanf("%d",&m);
printf("Cho chu ki T = ");
scanf("%f",&t1);
printf("So diem da cho n = ");
scanf("%d",&n);
for (i=1;i<=n;i++)
{
printf("x[%d] = ",i);
scanf("%f",&x[i]);
printf("y[%d] = ",i);
scanf("%f",&y[i]);
}
x[0]=1.0;
printf("%4cBANG SO LIEU\n",' ');
printf("%8cx%30cy\n",' ',' ');
for (i=1;i<=n;i++)
199
Tải về để xem bản đầy đủ
Bạn đang xem 20 trang mẫu của tài liệu "Giáo trình C++ - Chương 11: Nội suy và xấp xỉ hàm", để tải tài liệu gốc về máy hãy click vào nút Download ở trên
File đính kèm:
- giao_trinh_c_chuong_11_noi_suy_va_xap_xi_ham.pdf