Str. Calea Mărășești, nr. 157, Bacău, 600115 [608895]
c ROMÂNIA
MINISTERUL EDUCAȚIEI NAȚIONALE
UNIVERSITATEA „VASILE ALECSANDRI” DIN BACĂU
Facultatea de Științe
Str. Calea Mărășești, nr. 157, Bacău, 600115
Tel. +40 -234-542411, tel./ fax +40-234-571012
www.ub.ro ; e-mail: [anonimizat]
PROGRAMUL DE STUDII INFORMATICĂ
APLICATĂ ÎN ȘTIINȚĂ ȘI TEHNOLOGIE
LUCRARE DE DISERTA ȚIE
Coordonator științific:
Prof. univ. dr. Mihai Talmaciu Absolvent: [anonimizat]
2017
c
ROMÂNIA
MINISTERUL EDUCAȚIEI NAȚIONALE
UNIVERSITATEA „VASILE ALECSANDRI” DIN
BACĂU
Facultatea de Științe
Str. Calea Mărășești, nr. 157, Bacău, 600115
Tel. +40 -234-542411, tel./ fax +40-234-571012
www.ub.ro ; e-mail: [anonimizat]
PROGRAMUL DE STUDII INFORMATICĂ
APLICATĂ ÎN ȘTIINȚĂ ȘI TEHNOLOGIE
APROXIMAREA FUNCȚIILOR ȘI
APLICAȚII
Coordonator științific:
Prof. univ. dr. Mihai Talmaciu Absolvent: [anonimizat]
2017
3
CUPRINS
CUPRINS ………………………….. ………………………….. ………………………….. ………………………….. … 3
ARGUMENT ………………………….. ………………………….. ………………………….. ……………………….. 4
1. INTRODUCERE ………………………….. ………………………….. ………………………….. …………… 5
2. APROXIMAREA PRIN INTERPOLARE ………………………….. ………………………….. ………… 6
2.1. Interpolarea polinomială Lagrange ………………………….. ………………………….. …………….. 6
2.2. Algoritmul Aitken ………………………….. ………………………….. ………………………….. ……….. 9
2.3. Evaluarea restului la interpolarea Lagrange ………………………….. ………………………….. . 15
2.4. Diferențe divizate ………………………….. ………………………….. ………………………….. … 16
2.5. Formula lui Newton de interpolare ………………………….. ………………………….. …………… 23
2.6. Diferențe finite ………………………….. ………………………….. ………………………….. ………….. 24
2.7. Formule de interpolare pe noduri echidistante ………………………….. ……………………….. 25
2.8. Interpolarea polinomială Hermite ………………………….. ………………………….. …………….. 27
2.9. I nterpolarea prin funcții spline ………………………….. ………………………….. …………………. 29
3. APROXIMAREA ÎN MEDIE PĂTRATICĂ ………………………….. ………………………….. . 32
3.1. Problema fundamentală a aproximării liniare ………………………….. …………………………. 32
3.2. Teoreme fundamentale ………………………….. ………………………….. ………………………….. . 33
3.3. Aproximarea discretă în sensul celor mai mici pătrate ………………………….. ………. 34
4. APLICAȚII ………………………….. ………………………….. ………………………….. ………………… 36
5. BIBLIOGRAFIE ………………………….. ………………………….. ………………………….. …………. 64
4
ARGUMENT
Această lucrare este structurată pe cinci capitole, și anume:
Introducere;
Aproximarea prin interpolare;
Aproximarea prin medie pătratică;
Aplicații;
Bibliografie
În lucrare sunt prezentate câteva metode de aproximare a funcțiilor prin polinoame, cum ar
fi:
metode de aproximare prin interpolare;
metode de aproximare prin medie pătratică.
În primul capitol intitulat Introducere am prezentat problema aproximării unei funcții cu
ajutorul unui polinom, care se justifică prin modul simplu computerizat de evaluare a valorii
polinomului într -un punct.
La sfârșitul capitolului am făcut o clasificare a metodelor de aproximare a funcțiil or prin
polinoame.
În capitolul II, denumit, Aproximarea prin interpolare am prezentat metodele de
aproximare prin interpolare, și anume:
interpolarea polinomială Lagrange;
interpolarea polinomială Aitken;
interpolarea polinomială Newton;
interpolarea polinomială Hermite;
interpolarea prin funcții spline.
În cadrul acestor metode s -au prezentat și algoritmii de calcul, iar la primele trei metode și
listingul programelor în C++.
În capitolul III , intitulat, Aproximarea pr in medie pătratică am prezentat: Problema
fundamentală a aproximării liniare, Teoreme fundamentale și Aproximarea discretă în sistemul
celor mai mici pătrate.
În capitolul IV, denumit Aplicații am prezentat două aplicații de interpolare a funcțiilor prin
polinoame utilizând metodele Lagra nge și Newton care cuprinde pe lângă metoda matematică
de rezolvare și algoritmul de rezolvare împreună cu programul în limbajul de programare C++.
În încheiere am enumerat bibliografia de unde m -am inspirat în vederea elaborării lucrării
de disertație.
5
1. INTRODUCERE
În acest capitol am utilizat următoarele referințe bibliografice: ( [2], [6], [9] ).
Pentru o funcție
R],[:baf problema aproximării ei printr -un polinom se pune fie când
este dificil de evaluat f, fie când nu se cunoaște expresia analitică a lui f ci doar valorile ei în
anumite puncte
],[ba xi ,
n i,0 , obținute în general ca urmarea a unor măsurări și prezentate
într-un tabel. Mulțimea de puncte
],[ba xi ,
n i,0 , cu proprietatea
b x x xan …1 0
, o vom nota prin
],[bad și o vom numi diviziune a intervalului
],[ba . Punctele
ix ,
n i,0 vor
fi numite nodurile diviziunii.
Alegerea unui polinom pentru aproximarea funcției f se justifică prin modul simplu
computerizat de evaluare a valorii polinomului într -un punct.
Weierstrass, în 1885, demonstrează că „orice funcție continuă pe un interval
],[ba este
limita uniformă pe
],[ba a unui șir de polinoame” și dă un procedeu pentru calculul acestor
polinoame. Bernstein, în 1912, dă o metodă prin care putem calcula, pentru o funcție continuă
pe
],[ba , un șir
nnP)( de polinoame uniform convergent către f.
Fiind dată o funcție arbitrară
R],[:baf se pune problema stabilirii unui criteriu de
alegere a polinomului P care să aproximeze funcția dată. Astfel este necesar un instrument
matematic pentru măsurarea distanței dintre două funcții, care impune situarea într -un spațiu
vectorial normat complet și în funcție de stabilirea criteriului de alegere a polinomului P care
să aproximeze funcția f vom a vea metode de aproximare a funcțiilor prin polinoame și anume:
aproximarea prin interpolare ;
aproximarea în medie pătratică.
6
2. APROXIMAREA PRIN INTERPOLARE
Fie
R],[:baf și diviziunea
],[bad
:
b x x xan …1 0 .
Presupunem cunoscute valorile funcției f în nodurile
ix ,
n i,0 ale diviziunii
],[bad și
(eventual) derivatele de anumite ordine ale funcției f în aceste noduri. Se pune problema
determinării unui polinom P care are aceleași valori ca și f în noduri și (eventual) derivatele de
anumite ordine ale lui f și P în noduri să coincidă. Un astfel de polinom P îl vom numi polinom
de interpolare atașat fun cției f și diviziunii
],[bad .
2.1. Interpolarea polinomială Lagrange
În acest subcapitol am utilizat următoarele referințe bibliografice: ([ 3], [9]).
Fie
R],[:baf și
],[ba xi ,
n i,0 astfel încât
ix să fie toate distincte. Presupunem
că se cunosc valorile funcției f în nodurile
ix ,
)(i i xf f ,
n i,0 . Se pune problema
determinării unui polinom P astfel încât
i i f xP)( pentru
n i …,,1,0 .
Teorema 1. Există un polinom și numai unul de grad cel mult n astfel ca
i i f xP)(
pentru
n i …,,1,0 .
Demonstrație . Dacă P și Q sunt două polinoame de grad cel mult n astfel ca
i i i f xQ xP )( )(
pentru
n i …,,1,0 atunci polinomul
QPR este de asemenea de grad
cel mult n și
0)(ixR pentru
n i …,,1,0 . Deoarece polinomul R se anulează de cel puțin ( n
+ 1) ori și este de grad cel mult n, urmează că R este polinomul identic nul, adică
QP .
Existența polinomului P o demonstrăm construind efectiv polinomul P răspunzând
condițiilor
i i f xP)( pentru
n i …,,1,0 . Să găsim un polinom
k în spațiul vectorial al
polinoamelor de grad cel mult n, notat
nP , astfel ca
0)(kix
pentru
ki ,
1)(kkx .
Acest polinom se anulează de n ori pentru
1 10 …,,,kx xx ,
n k x x …,,1 . Deci
) )…( )( )…( )( ()(1 1 1 0 n k k k xx xx xx xxxx x
.
Deoarece
1)(kkx urmează că
n
kjj j kj
kx xxx
x
0)(
.
7
Cele ( n + 1) polinoame
k (
n k …,,1,0 ) sunt liniar independente în spațiul vectorial
al polinoamelor de grad cel mult n. Într -adevăr, fie o combinație liniară nulă:
0)( …)( )(11 00 x x xnn
.
Pentru
kxx ,
n k …,,1,0 , avem
0)( )( …)( )(11 00 k nn k kk k k x x x x
.
Deci
0k pentru
n k …,,1,0 .
Cele ( n + 1) polinoame
k formează deci o bază a acestui spațiu vectorial de dimensiune
n + 1. Polinomul P căutat se obține ca o combinație liniară de polinoame
k , adică polinomul
dat de formula Lagrange:
n
kkk x f xP
0)( )(
.
Într-adevăr,
in
kikk i f x f xP
0)( )(
, pentru
n i …,,1,0 .
Pentru calculul valorii polinomului de interpolare Lagrange într -un punct a folosind
relația
n
kkk x f xP
0)( )(
Algori tmul 1.
Intrări : n = gradul polinomului de interpolare
a = valoarea argumentului polinomului
x = tabloul absciselor punctelor de interpolare
y = tabloul ordonatelor punctelor de interpolare
Ieșiri:
= valoarea polinomului de interpolare în a
{
0
pentru
n i :0
{
1P
pentru
n j :0
dacă
ij atunci
) /() ( j i j xx xaP P
8
Pyi }
}.
Complexitatea metodei este
)(2nO , mai exact
24n operații (
)3 2(nn adăugări,
)12)(1( n n
înmulțiri și
)1(n împărțiri).
În continuare dăm o altă manieră de scriere a polinomului Lagrange. Fie w polinomul
de grad
)1(n definit prin:
) )…( )( ()(1 0 nxx xxxx xw
.
Atunci
) )…( )( )…( )( ()('1 1 1 0 n k k k k k k k k x x x x x x x xx x xw
.
Deci
)(') ()()(
k kkxwxxxwx
cu
1)(kkx .
Înlocuind
k prin expresia sa în formula Lagrange, obținem:
n
k k kkn
k k kk
xwxxfxwxwxxxwfxP
0 0 )(') ()()(') ()()(
.
În cazul când
1kf pentru orice k avem relația
1)(
0
n
kkx . Cum
1)(') (1)( )(
0 0
n
k k kn
kkxwxxxwx
obținem formulele baricentrice
n
k k kn
k k kk
xwxxxwxxf
xP
00
)(') (1)(') ()(
.
În acest caz numărul de operații este de ordinul
22n , dar o nouă evaluare a lui
)(xP ,
dacă se conservă valorile lui
)('1
kxw , nu va cere decât 2 n operații.
9
2.2. Algoritmul Aitken
În acest capitol am utilizat următoarele referințe bibliografice: ([6], [9]).
Atunci când se dorește modificarea punctelor de interpolare, se preferă următorul algoritm.
Lema 2.2.1. Fie Q și R două polinoame de interpolare de grad cel mult n, construite
respectiv pe
}, ,,{1 10 a …, xxxn și pe
}, ,,{1 10 b …, xxxn astfel încât:
i i i f xR xQ )( )( pentru
1 …,,1,0 n i
)(aQ
)(bR .
Atunci polinomul de interpolare de grad n + 1, constr uit pe
},, ,,{1 10 ba …, xxxn astfel
încât
i i f xP)( pentru
1 …,,1,0 n i
)(aP
)(bP
poate să se obțină plecând de la Q și R astfel:
abxRxa xQxbxP)() ()() ()(
.
Demonstrație. P este un polinom de grad cel mult n + 1.
iii ii
i fabfxa fxbxP ) ( ) ()(
,
1 …,,1,0 n i
ababaP) ()(
abbabP) ()(
.
În continuare se va construi polinomul P de grad cel mult n astfel încât:
i i f xP)(
pentru
n i …,,1,0 .
Considerăm următorul tabel triunghiular:
0 2 j j + 1 n
0
0,0P
x x0
1
0,1P
xx1
2
0,2P
2,2P
x x2
3
0,3P
2,3P
x x3
k
0,kP
2,kP
jkP,
1,jkP
x xk
n
0,nP
2,nP
nnP,
x xn
10
Fie următoarea relație de recurență:
pentru
n k …,,1,0
k k f P0, ;
pentru
1 …,,1,0 n j
pentru
n j jk …,,2,1
j kjk j jj k
jkx xx Px x x Px x
x P
)( ) ()( ) (
)(, ,
1,
.
Teorema 2.2.2. Polinomul de interpolare este egal cu
nnP, .
Demonstrație . Se va arăta că
jkP, cu
jk este polinomul de interpolare pe punctele
}, ,,{ 1 10 k jx …, xxx
, ceea ce va implica faptul că
nnP, este polinomul căutat.
Facem demonstrația prin recurență pe al doilea indice
k k f P0,
pentru
n k …,,1,0 .
Dacă ipoteza este verificată pentru j fixat atunci pentru
n jjk …,,1, :
jkP, este polinomul de interpolare pe punctele
}, ,,{ 1 10 k jx …, xxx ;
jjP, este polinomul de interpolare pe punctele
}, ,,{ 1 10 j jx …, xxx .
Aplicând lema 2.2.1.,
1,jkP este polinomul de interpolare pe punctele
},, ,,{ 1 10 kj j xx …, xxx
.
Listingul algoritmului Neville -Aitken
#include <stdio .h>
#include <conio .h>
#include <math .h>
#include <stdlib .h>
#include <iostream >
using namespace std;
#define eps 0.001
#define pi 3.14
float function(float x)
{
float y;
y=exp(x)*sqrt(x*x+1);
return y;
}
11
int main()
{
float z,x[100],y[100],t[100],a=0.5,b=4.0,p,t1,z1,er,d;
int n,i,k,j,m=1;
printf(" \nIntrodu punctul z in care se calculeaza valoarea polinomului de
interpolare(valoare de pe intervalul[0.5,4]): \n");
scanf("%f",&z);
printf("%d)",m);
printf("Introdu gradul polinomului de interpolare: \n");
scanf("%d",&n);
printf(" \n");
for(i=0;i<=n;i++)
{
t[i]=-cos((2*i+1)*pi/(2*n+2));
x[i]=(b+a)/2+t[i]*(b -a)/2;
}
for(i=0;i<=n;i++)
{
y[i]=function(x[i]);
}
for(i=0;i<=n;i++)
{
printf("x%d=%f f(x%d)=%f \n",i,x[i],i,y[i]);
}
z1=function(z);
er=-1; t1=y[0];
for(i=1;i< =n;i++)
{
for(k=i;k>=1;k –)
{
j=k-1;
y[j]=y[k]+(y[k] -y[j])*(z -x[i])/(x[i] -x[j]);
}
d=fabs(y[0] -t1);
t1=y[0];
if(d<eps){
er=eps;
p=y[0];
break;
}
else{
12
if((er<0)||(d<er)){
er=d;
p=y[0 ];
}
}
}
printf(" \n\nPunctul in care se calculeaza valoarea polinomului de interpolare: z=%f \n",z);
printf("Valoarea polinomului de interpolare in punctul z este: %f \n",p);
printf("Valoarea functiei in punctul z este: f(z)=%f \n",z1);
printf("Cea mai buna exactitate posibila: %f \n",er);
double z2,p1,x1[100],y1[100],pas=0,er1,d1,t2;
int n1,k1;
m=1;
printf(" \n\n\n%d)",m);
printf("Introdu gradul polinomului de interpola re:\n");
scanf("%d",&n1);
printf(" \n");
pas=(b -a)/n1;
for(i=0;i<=n1;i++)
{
x1[i]=a+(pas*i);
}
for(i=0;i<=n1;i++)
{
y1[i]=function(x1[i]);
}
for(i=0;i<=n1;i++)
{
printf("x%d=%f f(x%d)=%f \n",i,x1[i],i,y1[i]);
}
z2=function(z);
er1=-1; t2=y1[0];
for(i=1;i<=n1;i++)
{
for(k1=i;k1>=1;k1 –)
{
j=k1-1;
y1[j]=y1[k1]+(y1[k1] -y1[j])*(z -x1[i])/(x1[i] -x1[j]);
}
13
d1=fabs(y1[0] -t2);
t2=y1[0];
if(d1<eps){
er1=eps;
p1=y1[0];
break;
}
else{
if((er 1<0)||(d1<er1)){
er1=d1;
p1=y1[0];
}
}
}
printf(" \n\nPunctul in care se calculeaza valoarea polinomului de interpolare: z=%f \n",z);
printf("Valoarea polinomului de interpolare in punctul z este: %f \n",p1);
printf("Valoarea functiei in punctul z este: f(z)=%f \n",z2);
printf("Cea mai buna exactitate posibila: %f \n",er1);
getch();
}
14
15
2.3. Evaluarea restului la interpolarea Lagrange
În acest subcapitol am utilizat următoarele referințe bibliografice: ([ 1], [6], [ 8], [9]).
Fie
R],[:baf și fie
],[ …,,,10 ba x xxn . Ce se poate spune despre restul
)()( xPxf
, unde P este polinomul de interpolare, astfel că
i i f xP)( ,
n i …,,1,0 .
Se remarcă faptul că x nu poate să aparțină lui
],[ba . În aceste caz, se spune că avem o
extrapolare.
Aproximarea realizată cu ajutorul polinomului de interpolare Lagrange este nelocală, în
sensul că la evaluarea funcției de interpolat contribuie informații din toate punctele de
interpolare, nu numai din cele din vecinătatea argumentului considerat. Aceasta conduce la
oscilații complet străine de compor tarea reală a funcției aproximate. Un exemplu în acest sens
este pentru aproximarea funcției
x xf /1)( . Distribuirea în mod echidistant a punctelor de
interpolare reduce oscilațiile până aproape la dispariție.
Teorema 2.3.1.
Dacă
],[1ba Cfn și
b x x xan …1 0 atunci
)()!1() )…( )( ()( )()1( 1 0n nfnxx xxxxxPxf
,
unde
b xx xx an ), max( ), min(0 .
Demonstrație. Fie
R],[:baF , definită prin
)( )()( )( ywcyPyf yFn
unde c este o constantă aleasă astfel încât
0)(xF .
Deci
)(/))()(( xwxPxf c , unde
) )…( )( ()(1 0 nxx xxxx xw .
Din condițiile de interpolare
0)(ixF ,
n i …,,1,0 .
Deci F are n + 2 zerouri distincte și anume:
nx xxx …,,,,10 pe
],[ba . Aplicând teorema
Rolle pe acest interval,
'F va avea cel puțin
1n zerouri distincte,
"F cel puțin
n zerouri
distincte. Din aproape în aproape,
)1(nF va avea cel puțin un zero
),(ba și deci
)!1()( )( 0)1( )1( nc f Fn n
.
Sub această formă restul nu poate fi determinat în practică, deoarece punctul
, care
depinde de x și de nodurile de interpolare nu este cunoscut. Dacă, însă, se cunoaște estimația
|)( | sup)1(
],[y f Mn
bayn
, se poate evalua
16
)!1(|)(||)( )(|nxw MxPxfn.
Interpolarea polinomială Lagrange este utilă pentru construirea formulelor de aproximare
a integralelor, de rivatelor și ecuațiilor diferențiale.
În practică, pentru a aproxima o mulțime de puncte de interpolare, se va limita la
polinoame de grad cel mult 10.
Se va prefera pentru polinoame de grad mai mare decât 10, să se utilizeze funcții spline
care au cel e mai bune proprietăți de stabilitate.
Exemplul 1. Să se construiască polinomul de interpolare Lagrange a funcției f definită
pe
]5,5[ prin
2
)(xexf .
Rezolvare. Se constată că se aproximează curba
2
)(xexf la centrul intervalului,
iar la extremități se produc oscilații numite efecte de margine.
Exemplul 2. Să se construiască polinomul de interpolare Lagrange a funcției f definite
pe
]5,5[ prin
211)(
xxf
.
Rezolvare. Se remarcă faptul că se produc efecte de margine.
Exemplul 3. Se dă funcția f prin următorul tabel:
ix
0 1 2
)(ixf
4 1 10
Se cere polinomul Lagrange
)(xP .
Rezolvare.
2( 1)( 2) ( 0)( 2) ( 0)( 1)( ) ( 4) 1 10(0 1)(0 2) (1 0)(1 2) (2 0)(2 1)
2 3 4.x x x x x xPx
xx
2.4. Diferențe divizate
În acest subcapitol am utilizat următoarele referințe bibliografice: ([ 1], [2], [3], [7]).
Formula restului amintește oarecum de restul unei alte aproximări polinomiale a unei
funcții de n ori continuu diferen țiabilă – formula lui Taylor. Vom arăta că și polinomul lui
Lagrange se poate exprima sub o formă similară cu polinomul lui Taylor, în care derivatele se
înlocuiesc cu diferențe divizate.
17
Diferența divizată de ordinul zero în
0x este
)( ][ 0 0 xf xf
.
Prin definiție diferența divizată de ordinul întâi pe nodurile distincte
ix și
jx este numărul
notat cu
fjixx ],[ și este dat de:
i ji j
fjix xxf xf
xx
)( )(
],[
.
Diferența divizată de o rdinul al doilea pe nodurile distincte
ix ,
jx și
kx este:
i kfji fkj
fkjix xxx xx
xxx
],[ ],[
],,[
.
În general, diferența divizată de ordinul k pe
1k noduri distincte se definește prin
recurență:
1 11 1 2
1 21] ,…,[ ] ,…,[
] ,…,,[x xx x x x
x xx
kfk f k
f k
în care diferențele divizate de la numărător sunt de ordin
1k .
Rezultatul care urmează precizează că diferența divizată depinde de nodurile în care este
definită și de valorile funcției f în aceste noduri.
Teorema 2.4.1. Pentru orice
2k și orice noduri distincte
kx xx ,…,,21 are loc formula:
12
1
1()[ , ,…, ]
()ki
kfki
ij
j
jifxx x x
xx
.
Demonstrație. Facem inducție după k. Pentru
2k , formula din enunț se verifică conform
definiției diferenței divizate de ordinul întâi. Presupunem adevărată formula din enunț pentru k
noduri distincte oarecare și verificăm pentru
1k noduri. Utilizând formula de recurență din
definiț ia diferenței divizate de ordinul k și ipoteza de inducție avem:
11
1 2 1121 111
2 1 1
11 1
112
1 1 1
12( ) ( ) ( ) 1[ , ,…, ]
( ) ( ) ( )
( )( ) ( )( ) ()
( ) ( ) ( )kki i k
kfk k kii ki j i j k j
j j j
j i j i
ki i i i k
kki
k i j j
jj
jif x f x f xx x xxxx x x x x x
f x x x f x x x fx
x x x x x x
1
1
1().
()ki
ki
ij
j
jifx
xx
18
Diferențele divizate definite prin recurență furnizează un algoritm de complexitate
)(2nO
.
Algoritmul 2.4.
Intrări : n + 1 = numărul de puncte
x = tabloul celor n + 1 puncte
y = tabloul valorilor funcție f în cele n + 1 puncte
Ieșiri: y = diferențele divizate de ordin 0 .. n în
0x
{
pentru
n i :1
pentru
inj :
) /() ( 1 1 j j j j j xx y y y
}.
Listingul p olinomul ui Newton cu diferențe divizate
#include <stdio .h>
#include <conio .h>
#include <math .h>
#include <stdlib .h>
#include <iostream >
using namespace std;
#define eps 0.001
#define pi 3.14
float function(float x)
{
float y;
y=exp(x)*sqrt(x*x+1);
return y;
}
int main()
{
float z,x[100],y[100],t[100],a=0.5,b=4.0,p,z1,E;
int n,i,k,m=1;
19
printf(" \nIntrodu punctul z in care se calculeaza valoarea polinomului de
interpolare(valoare de pe intervalul[0.5,4]): \n");
scanf("%f",&z);
do{
printf("%d)",m);
printf("Introdu gradul polinomului de interpolare: \n");
scanf("%d",&n);
printf(" \n");
for(i=0;i<=n;i++)
{
t[i]=-cos((2*i+1)*pi/(2*n+2));
x[i]=(b+a)/2+t[i]*(b -a)/2;
}
for(i=0;i<=n;i++)
{
y[i]=function(x[i]);
}
for(i=0;i<=n;i++)
{
printf("x%d=%f f(x%d)=%f \n",i,x[i],i,y[i]);
}
z1=function(z);
for(i=1;i<=n;i++)
{
for(k=n;k>=i;k –)
{
y[k]=(y[k] -y[k-1])/(x[k] -x[k-i]);
}
}
p=y[n];
for(i=n -1;i>=0;i –)
{
p=y[i]+(z -x[i])*p;
}
E=z1 -p;
if(fabs(E)>eps){ printf(" \nEroarea de calcul in punctul z a polinomului de interpolare este
mai mare decit precizia data epsilon=%f \n",eps);
20
printf(" \nIntrodu un grad mai mare pentru polinomul de interpolare \n");
}
m++;
}while(fabs(E)>eps);
printf(" \nPunctul in care se calculeaza valoarea polinomului de interpolare: z=%f",z);
printf(" \nValoarea polinomului de interpolare in punctul z este: %f",p);
printf(" \nValoarea functiei in punctul z este: f(z)=%f",z1);
printf(" \nEroarea in punctul z este: %f",fabs(E));
float z2,p1,E1,x1[100],y1[100],pas=0;
int n1;
m=1;
do{
printf(" \n\n%d)",m);
printf("Introdu gradul polinomului de interpolare: \n");
scanf("%d",&n1);
printf(" \n");
pas=(b -a)/n1;
for(i=0;i<=n1;i++)
{
x1[i]=a+(pas*i);
}
for(i=0;i<=n1;i++)
{
y1[i]=function(x1[i]);
}
for(i=0;i<=n1;i++)
{
printf("x%d=%f f(x%d)=%f \n",i,x1[i],i,y1[i]);
}
z2=function(z);
for(i=1;i<=n1;i++)
{
for(k=n1;k>=i;k –)
{
y1[k]=(y1[k] -y1[k-1])/(x1[k] -x1[k-i]);
}
}
21
p1=y1[n1];
for(i=n1 -1;i>=0;i –)
{
p1=y1[i]+(z -x1[i])*p1;
}
E1=z2 -p1;
if(fabs(E1)>eps){ printf(" \nEroarea de calcul in punctul z a polinomului de interpolare
este mai mare decit precizia data epsilon= %f\n",eps);
printf(" \nIntrodu un grad mai mare pentru polinomul de
interpolare \n");}
m++;
}while(fabs(E1)>eps);
printf(" \nPunctul in care se calculeaza valoarea polinomului de interpolare: z=%f",z);
printf(" \nValoarea polinomului de interpolare in punctul z este: %f",p1);
printf(" \nValoarea functiei in punctul z este: f(z)=%f",z2);
printf(" \nEroarea in punctul z este: %f",fabs(E1));
getch();
}
22
23
2.5. Formula lui Newton de interpolare
În acest subcapitol am utilizat următoarele referințe bibliografice: ([ 2], [3], [8], [9]).
Cu ajutorul diferenței divizate putem exprima restul la interpolarea Lagrange sub o altă
formă. Avem:
,)
) () ()(
) ()(() ()( )( )( )(
0
000 0
n
i
ijj i ii
n
iin
iin
in
ijj j ij
i
xx xxxf
xxxfxxxxxx
xf xf xPxf
de unde, conform relației din teorema 5.2.4.1.
fn n x xxxwxPxf ] ,…,,[)( )()( 0
,
unde (am renotat)
) () )( ()(1 0 n n xx xxxx xw .
Polinomul Lagrange
)( )( xLxP n
se scrie:
))( )((…))( )(()( )(1 1 2 1 x LxL xLxL xLxLnn n
.
Pentru fiecare m,
nm2 , polinomul
)( )(1x LxLm m
are gradul cel mult m și
0)( )( 1 j m j m x L xL
,
pentru
mj0 .
Deci
)( )( )(1 1 1 x w ax LxLm m m m
,
unde
1
01 ) ( )(m
jj m xx x w
.
Pentru
mxx obținem
)( )( )(1 1 1 m m m m m m x w a x L xf ; comparând cu expresia
lui
)()( xPxf , pentru
1mn și
mxx , deducem
fm f m m m x xx x xxx a ] ,…,,[ ] ,…,,,[ 10 1 10 1
.
Se obține polinomul lui Newton de interpolare:
.) )…( )( (] ,…,,[…) (],[)( )( 1 1 0 10 0 10 0 n fn f xx xxxx x xx xx xx xf xP
24
Calculul polinomului Newton de interpolare se exprimă prin:
Algoritmul 2.5.
Intrări : n = gradul polinomului de interpolare
a = abscisa în care se calculează polinomul
x = tabloul punctelor de interpolare
y = tabloul valorilor prin f a punctelor de interpolare
Ieșiri:
= valoarea po linomului de interpolare în a
{
apelează Dif._div. (n,x,y)
0y s
1P
pentru
n i :1
{
) (1ixaP P
iyPss
}
s }.
2.6. Diferențe finite
În acest subcapitol am utilizat următoarele referințe bibliografice: ([ 6], [8], [10], [12]).
Nodurile de interpolare
nx xx …,,,10 sunt echidistante dacă
ihx xi0 ,
n i …,,1,0 .
Numim diferență finită de ordinul întâi în x
)( ) ( )( xf hxf xf .
Se verifică direct că
este un operator liniar, iar diferența finită de ordinul k în x se
definește prin relația
))( ( )(1xf xfk k
,
…,3,2 k .
Formula de calcul pentru
)(xfk este următoarea:
) ( )1( )(
0jhxfC xfk
jj
kjk k
.
Pentru demonstrație, procedăm prin inducție. Avem:
). ( )1( ) ( )1()( ) ( )(
0 01
jhxfC jhhxfCxf hxf xf
k
jj
kjkk
jj
kjkk k k
Făcând în prima sumă schimbarea de indice
1 'jj și renotând apoi
jj' , obținem:
25
). ( )1()0( )1() () ( )1())1(( )1() ( )1( )' ( )1( )(
1
0110
101
11 11
1)1(1011
1'1' 1' 1
jhxf Ch xf C jhxfC Ch kxf CjhxfC hjxf C xf
k
jj
kj kkkk
jj
kj
kj kk
kk kk
jj
kjkk
jj
kjk k
Prin
inducție după k, se arată că are loc:
kik
fki i ihkxfx xx!)(] ,…,,[1
.
2.7. Formule de interpolare pe noduri echidistante
În acest subcapitol am utilizat următoarele referințe bibliografice: ([ 2], [6], [ 7], [9]).
Să presupunem că punctul x în care se face interpolarea pe nodurile echidistante
nx xx …,,,10
este apropiat de
0x . Făcând transformarea
thxx0 în polinomul lui Newton
de interpolare și ținând cont de legătura dintre diferențe finite și divizate, obținem după
simplificări:
2
0 0 0 0 0( 1) ( 1)…( 1)( ) ( ) ( ) ( ) … ( ).2! !n t t t t t nP x th f x t f x f x f xn
Aceasta se numește formula lui Newton progresivă pe noduri echidistante.
Dacă x este apropiat de
nx , făcând transformarea
i ni 1 ' ,
n i …,,1,0 , formula lui
Newton devine:
). )…( )( (] ,…, ,[…) )( (] , ,[) (] ,[)( )(
1 1 0 11 2 1 1
xx xxxx x xxxxxx x xx xx xx xf xP
n n f nnn n f n nn n f nn n
Făcând transformarea
th xxn și ținând cont de legătura dintre diferențe finite și
divizate, obținem după simplificări:
2
12
0( 1)( ) ( ) ( ) ( ) …2!
( 1)…( 1)( ),!n n n n
nttP x th f x t f x f x
t t t nfxn
numit polinomul Newton regresiv pe noduri echidistante.
Formula lui Newton progresivă pe noduri echidistante conduce la următorul algoritm.
26
Algoritmul 2.7.
Intrări : n + 1 = numărul de puncte
h = pasul
x = tabloul punctelor de interpolare
ah x0 = abscisa în care se calculează polinomul
Ieșiri:
= valoarea polinomului de interpolare în
ah x0
{
)(0xf
pentru
n k :1
{
1P
pentru
1:1 k q
q qaP P /)1 (
dacă k este par atunci
)(0xf d
altfel
)(0xf d
pentru
k j :1
{
1c
pentru
j q :1
q qkcc /)1 (
dacă
) ( jk este par atunci
) (0jh xfcd d
altfel
) (0jh xfcd d
}
}
}.
27
Exemplu. Se dă funcția f prin următorul tabel:
ix
0 1 2 3
iy
4 1 10 29
Se cere polinomul Newton progresiv.
Rezolvare. Se întocmește un tabel cu diferențe finite și se obține:
)2 )(1(!36)1(!2454 )1 0( t tt tt t t P
45 )(2 3 t ttxP
.
2.8. Interpolarea polinomială Hermite
În acest subcapitol am utilizat următoarele referințe bibliografice: ([6], [ 8], [9] )
Fie
R],[:baf și
nx xx …,,,10 noduri de interpolare, toate distincte. Se caută un
polinom P astfel încât:
i i ii i if xf xPf xf xP
' )(' )(')( )(
,
n i …,,1,0 .
Teorema 2.8.1. Există un polinom P și numai unul singur de grad cel mult
12n astfel
încât:
i ii if xPf xP
' )(')(
pentru
n i …,,1,0 .
Demonstrație . Presupunem că există polinoamele P și Q de grad cel mult
12n astfel
ca
i i i f xQ xP )( )( și
i i i f xQ xP ' )(' )(' pentru
n i …,,1,0 . Polinomul
QPR , care
are de asemenea gradul cel mult
12n , admite fiecare
ix ca rădăcină dublă (
0)(' )( i i xR xR
). Deoarece R are cel puțin
2 2n rădăcini și are gradul cel mult
12n
urmează că R este polinomul nul și deci
QP .
Să arătăm că polinomul
n
ii in
ii i fxk fxh xP
0 0')( )( )(
,
unde:
)()](') (21[)(2x x xx xhiiii i
)() ()(2x xx xkii i
cu
28
n
ijj j ij
iixxxx
x
0)(
îndeplinește condițiile din enunț.
Deoarece
i este un polinom de grad n,
ih și
ik sunt polinoame de grad
12n . Cum
ij jix)(
,
1)()](') (21[)(2 iiiii i ii x x xx xh
,
0)() ()(2 iii i ii x xx xk
,
pentru
ij avem:
0)()](') (21[)(2 jiiii j ji x x x x xh
,
0)() ()(2 jii j ji x x x xk ,
adică
k k f xP)( .
)(')(2)](') (21[)()('2 )('2x x x xx x x xhi i iii iii i
,
)(')() (2)( )('2x x xx x xki ii i i
.
Deci:
0)('2)('2)(' ii ii ii x x xh
,
1)( )('2ii ii x xk ,
pentru
ij avem
0)('jixh ,
0)('jixk , deci
k k f xP ' )(' .
Pentru evaluarea erorii la interpolarea Hermite avem
Teorema 2.8.2.
Dacă
],[22ba Cfn și
b x x xan …1 0 atunci
)()!2 2() …() () ()( )()22(2 2
12
0 n nfnxx xx xxxPxf
,
unde
b xx xx an ), max( ), min(0 .
Demonstrația este asemănătoare celei de la interpolare Lagrange.
Polinomul P având expresia
n
ii in
ii i fxk fxh xP
0 0')( )( )(
se numește polinomul lui Hermite construit pe baza condițiilor de interpolare
i i ii i if xf xPf xf xP
' )(' )(')( )(
,
n i …,,1,0 .
29
Pentru calculul valorii polinomului Hermite într -un punct a a fost definit un
)(2nO –
algoritm.
Algoritmul 2.8.1.
Intrări : n + 1 = numărul nodurilor de interpolare
a = abscisa în care se calculează valoarea polinomului
x = tabloul celor n + 1 puncte de interpolare
y = tabloul valorilor funcției de interpolat în cele n+1 puncte
dy_ = tabloul valorilor derivatelor funcției de interpolat în cele n+1 puncte
Ieșiri:
h = valoarea polinomului Hermite în a
{
s
0
pentru
n i :0
{
1
0 _d
pentru
n j :0
dacă
ij atunci
{
) /() ( j i j xx xa
) /(1 _ _j ixx d d
s
s + +
i i i i dy xa y dlxa _ ) ( _) (21
}
h s
}.
2.9. Interpolarea prin funcții spline
În acest subcapitol am utilizat următoarele referințe bibliografice: ([ 2], [6], [9]).
Studiul interpolării polinomiale ne conduce la două constatări. Prima este costul calculelor
care devine ridicat când gradul polinomului este mare, iar a doua este că se „stăpânește” rău
efectul de margine dacă intervalul de interpolare este mare.
O altă aproximare este cea locală, adică făcând interpolarea pe subintervale prin funcții
spline.
30
Fie
nx xx …,,,10 punctele de interpolare. Pe fiecare subinterval
] ,[1iixx se caută un
polinom de interpolare
iS astfel încât condițiile de interpolare să fie îndeplinite:
i i ii f xf xS )( )(
pentru
n i …,,1,0 .
În plus se cere ca
iS să fie de clasă
2C . De asemenea, pentru
1 …,,2,1 n i :
)(' )('1 ii i i xS x S
, continuitatea primei derivate,
)(" )("1 ii i i xS x S
, continuitatea derivatei a doua.
Condițiile de interpolare și cele de clasă impun a căuta
iS de grad ul trei. Exprimăm
iS
în funcție de
if" și
1"if (necunoscute încă).
Prin interpolare liniară avem, punând
i i i x x h1 :
ii
i
ii
i ihxxfhx xf xS
11" " )("
.
Integrând de două ori obținem:
) ( ) (6) ("6) (" )( 13
13
1
i i ii
ii
i
ii
i i xxbx xahxxfhx xf xS
,
unde
ia și
ib sunt constante de integrare, care se pot calcula ținând cont de
i ii f xS)( și
1 1) (i ii f xS
, adică
i iii
i f hahf 6"2
și
12
16" i iii
i f hbhf ,
de unde
6"i
i
ii
ihfhfa
și
6"11 i
i
ii
ihfhfb .
Înlocuind
ia și
ib prin expresiile lor obținem:
). ( ) () (6 6) (" ) (6 6) (" )(
1
13
1 13
1
i
ii
i
iiii
ii
i ii
ii
i i
xxhfx xhfxxh
hxxf x xh
hx xf xS
Punem condiția ca derivata să fie o funcție continuă pe intervalul dat. Avem
ii
ii i
ii
ii
ii
i ihf
hf h
hxxfh
hx xf xS12
12
1
6 2) ("6 2) (" )('
)(' )('1i i ii x S xS ,
echivalent cu:
31
1 11 1 1
11
13"6"6"3"
ii
ii i
ii
i
ii
ii i
ii
ihf
hf hfhfhf
hf hfhf,
sau încă pentru
1 …,,2,1 n i :
11 1
1 1 1 1 6 " ") (2"
ii i
ii i
i i i i i i ihff
hf ffhfhh fh
.
Am obținut un sistem de
1n ecuații liniare cu
1n necunoscute. Precizând valorile
lui
0"f și
nf" se poate rezolva sistemul liniar tridiagonal.
Exemplu. Pe intervalul
]5,5[ fie funcțiile f și g definite prin
2
)(xexf și
11)(2
xxg
. Să se facă interpolarea prin funcții spline cu 3, 5 și respectiv 9 puncte
echidistante și să se dea valorile derivatelor la extremități.
Rezolvare. Se constată că pentru 9 puncte de interpolare avem o foarte bună
aproximare și nici un efect de margine.
32
3. APROXIMAREA ÎN MEDIE PĂTRATICĂ
3.1. Problema fundamentală a aproximării liniare
În acest subcapitol am utilizat următoarele referințe bibliografice: ([ 2], [6], [9]).
Exemplu. Să aproximăm pe
)(xf pe intervalul
],[ba printr -un polinom
)(xP .
Se pune problema să se găsească P astfel ca:
a)
1
02))( )(( dxxPxf să fie minimă;
b)
1
021
02))( )(( ))( )(( dxxPxfdxddxxPxf să fie minimă;
c)
|)( )(| max
1 0xPxf
x
să fie minimă.
Este necesară o măsură numerică a „distanței” dintre funcția de aproximat și polinomul
aproximant într -un spațiu vectorial normat.
Fie X un spațiu vectorial normat peste corpu l K (adesea
RK ). Notăm prin
),(
produsul scalar și prin
norma.
Fie
nx xx ,…,,21 , n elemente din X liniar independente. Trebuie să aproximăm y printr –
o combinație liniară de
nx xx ,…,,21 , adică, să găsim
na aa ,…,,21 cu
K ai , astfel ca
) … ( 22 11 nnxa xaxay
să fie cât mai mic posibil.
Definiția 3.1. Cea mai bună aproximare a lui y printr -o combinație liniară de
n iix 1)(
este un element
nnxa xaxa …22 11 astfel ca:
Kb xb xbxby xa xaxay i nn nn ,) … ( ) … ( 22 11 22 11
.
Cea mai bună aproximare constă în a găsi minimul normei erorii.
Drept norme se consideră:
1) dacă
RX atunci
x x (modulul);
2) dacă
nXR atunci
2/12 2
22
1 ) … (nx x x x (norma euclidiană);
3) dacă
),(baCX (= spațiul vectorial al funcțiilor continue pe
],[ba ) atunci:
a)
)( max xf f
bxa (norma convergenței uniforme);
b)
2/1
2|)(|
b
adx xf f (norma convergenței pătratice).
33
3.2. Teoreme fundamentale
În acest subcapitol am utilizat următoarele referințe bibliografice: ([ 2], [9]).
Teorema 3.2. Fie n elemente liniar independente din
nx xxX …,,,:21 . Pentru orice
Xy
, problema celei mai bune aproximări are o soluție.
Verificare. Trebuie arătat că funcția norma erorii
) … ( ) ,…,,( 22 11 21 nn n xa xaxay a aae
are un minim. Pentru aceasta, se arată că această funcție este continuă și că este suficient a
limita mulțimea de variație a lui
na aa ,…,,21 la o mulțime mărginită și închisă (exteriorul
acestei mulțimi nu conține minimul) pentru a aplica teorema: orice funcție continuă pe o
mulțime închisă și mărginită din
nR își atinge minimul său.
Corolarul 3.2.1. Fie
),(baCf , n întreg fixat.
Problema: să se găsească
)(ia astfel ca:
|) … ()(| max min 1 0n
nbxaaxa xa a xf
i
are o soluție. Această soluție este unică.
Aceasta se numește aproximarea Tchebycheff de grad n a lui f.
Corolar 3.2.2. Fie
),(baCf , n întreg fixat.
Problema: să se găsească
)(ia astfel ca:
b
an
nadx xa xa a xf
i|) … ()(| min1 0
are o soluție. Această soluție este unică.
Aceasta se numește aproximarea continuă în sensul celor mai mici pătrate.
Corolarul 3.2.3. Date fiind
kx xx ,…,,10 ,
1k puncte distincte cu
nk , problema: să se
găsească
)(ia astfel ca:
n
in
in i i iaxa xaxa a xf
i 02 2
2 1 0 )) … ()(( min
are o soluție.
Aceasta se numește aproximarea discretă în sensul celor mai mici pătrate.
34
3.3. Aproximarea discretă în sensul celor mai mici pătrate
În acest subcapitol am utilizat următoarele referințe bibliografice: ([ 2], [9]).
a) Cazul general
Să se determine
R R:P unde
)( …)( )( )(11 00 xua xuaxuaxPnn
,
ceea ce revine la a determine
)1(n numere
)(ia .
Cele
)1(n funcții liniar independente
)(iu sunt date. Se pot alege,
k
k xxu)( (
n k …,,1,0
) sau funcțiile trigonometrice sau cele exponențiale. Se cunosc
)1(k perechi
),( iiyx
(
k i …,,1,0 ) cu
nk , unde
ix sunt distincte. Fie
)(i i i xPyr ,
k i …,,1,0 .
Să se determine
1() Rn
i aa astfel ca
k
iir aR
02)( să fie minim.
Considerând norma euclidiană în
1kR , avem
||||)(2r aR , unde
1
0 R )(
k
ki irr .
Luăm:
1
0 R )(
k
ki iy y .
))(()( ji ij xu u U , matricea rectangulară cu
)1(k linii și
)1(n
coloane.
Avem
) , (),(||||2UayUay rr r .
Minimul lui R este caracterizat prin
) ( )( baRaR ,
}0{ R ,
1nbR , adică:
).,( ), (2 0) , () , ()) ( ), ( () , (
2UbUb UbUayUb UayUb Uay UayUaybaUybaUy UayUay
Pentru
0 , împărțind prin
2 și făcând
să tindă spre zero obținem
), ( 0 UbUay
. pentru
0 , împărțind prin
2 și făcând
să tindă spre zero obținem
), ( 0 UbUay
. Cele două inegalități dau
0), ( UbUay .
Notând cu
TU transpusa lui U obținem
0), ( bUaUyUT T ,
1nbR . Deci
valoarea minimă a funcției R verifică sistemul liniar
yU UaUT T , unde
UUT este o matrice
pătratică de ordin
1n .
Acest sistem se numește sistemul ecuațiilor normale.
35
b) Cazul liniar
În acest caz particular, se caută polinomul de grad unu,
xa axP1 0)( astfel ca
k
ii ik
ii i bxby axay
02
0 1
02
0 1 ) ( ) (
,
2
10),( R bb . Fie
1)(k
ix X R și
1)(k
iy Y R
. Sistemul de ecuații normale devine:
kk
kk
yyyy
x xx x aa
xxxx
x xx x
210
2 1 0 10
210
2 1 01 1 1 1
1111
1 1 1 1
,
adică
k
iiik
ii
k
iik
iik
ii
yxy
aa
x xx k
00
10
02
001 . Notând
k
iixkx
011 ,
k
iiyky
011 ,
deoarece
yxxyyxyxyyxxi i ii i i ) )( ( , avem
yxxkyykxyxkyyxxkk
iik
iik
iiik
ii i 0 0 0 0 1 1 11) )( (11
, care, în statistică, se
numește covarianța vectorilor X și Y și se no tează
k
ii i yyxxkYX
0) )( (11),(
.
Dacă
YX , obținem
k
iixxkX
02 2) (11)( .
Rezolvarea sistemului dă
)(),(
21XYXa
,
xay a1 0 .
36
4. APLICAȚII
În acest capitol am utilizat următoarele referințe bibliografice: ([6], [9], [10], [11], [12]).
4.1. Problema nr. 1
Să se găsească polinomul de interpolare Lagrange pentru funcția
3)( xxf și nodurile
5 ,3 ,13 2 1 x x x
.
Rezolvare.
. )( )( )()( )(
2 32
1 313
3 23
1 212
3 13
2 1211 1
x xxx
x xxxxfx xxx
x xxxxfxxxx
xxxxxfxxxx
xf xLn
in
ijj j ij
i n
Astfel,
2438)3 )(1(274)5 )(1(18)5 )(3()(3 x x x x x xxL
,
de unde
423782495)(2
3 x x xL
.
Algoritmul
Algoritmul pentru determinarea valorii polinomului Lagrange este
Intrări: n = numărul nodurilor de interpolare
a = punctul în care se calculează valoarea polinomului
x = vectorul cu nodurile de interpolare
f = funcția din metodă
Ieșiri:
= valoarea polinomului Lagrange în a
{
0
pentru i
1 : n
{
p
1
pentru j
1 : n
dacă
ij atunci
) /() ( j i j xx xap p
)(ixfp
}}
37
Programul sursă C ++
// Polinomul de interpolare Lagrange
#include <stdio .h>
#include <conio .h>
#include <iostream >
using namespace std;
#define max 10
float funct(float t)
{
float f;
f=t*t+t+1;
return f;
}
int main ( )
{
float l,p,a;
float x[max+1];
int n,i,j;
cout << "dati valoarea in care se doreste aproximarea functiei" << endl;
cin >> a; cout << "dati numarul de noduri" << endl;
cin >> n; cout << "dati nodurile" << endl;
for (i=1;i<=n;i++)
{
cout << "x[" << i << "]=";
cin >> x[i];
}
cout << endl;
l=0;
for (i=1;i<=n;i++)
{
p=1;
for (j=1;j<=n;j++)
if (j !=i)
p=p*(a -x[j])/(x[i] -x[j]);
l=l+funct(x[i])*p;
38
}
cout << "Val. pol. in " << a << " este " << l << endl;
getche();
}
39
4.2. Problema nr. 2
Să se scrie polinomul de interpolare Newton pentru funcția
1 )(2 x xxf și nodurile
de interpolare 1, 2, 5.
Polinomul de interpolare Newton este
3
21
121 1 ) ( ] ,…,,[ )( )(
kk
fk xx x xx xf xP
.
Astfel
)2 )(1 )(821
13
21()1 )(13
21(1) )( ()
) ()(
) ()(
) ()(() ()
) ()(
) ()((1) )( (
) ()() (
) ()(1) ( ],,[) ( ],[)( )(
2 12
3133
2
2122
2
111112
2122
2
11112 13
13
112
12
12
13211
121 1
x x xxxxx
x xxf
x xxf
xxxfxx
x xxf
xxxfxxxx
xxxfxx
xxxfxx xxx xx xx xf xP
jjj
jjj
jjjjjj
jjji
ijjj ii
i
ijjj iif f
de unde rezultă
)2 )(1(841)1(251)( x x x xP
.
Algoritmul
Algoritmul asociat metodei Newton este
Intrări: n = numărul nodurilor de interpolare
a = numărul în care se calculează valoarea polinomului
x = vectorul cu nodurile de interpolare
f = funcția din metodă
Ieșiri:
1 = valoarea polinomului Newton în a
{
2s
0
pentru k
2 : n
{
p1
1
pentru
1 : k – 1
)(ixfptt
40
s1
0
pentru i
1 : k
{
p2
1
pentru j
1 : k
dacă
ij atunci
) (2 2 j ixx p p
2/)( 1 1 pxfs si
}
212 2 pss s
}
2)( 11s xf
}
Programul sursă C ++
// Polinomul de interpolare Newton
#include <stdio .h>
#include <conio .h>
#include <iostream >
using namespace std;
#define max 10
float funct(float t)
{
float f;
f=t*t+t+1;
return f;
}
int main ( )
{
float l,s,p,p1,a;
float x[max+1];
int n,i,j,k,m;
cout << "dati numarul de noduri" << endl;
41
cin >> n;
cout << "dati nodurile" << endl;
for (i=1;i<=n;i++)
{
cout << "x[" << i << "]=";
cin >> x[i];
}
cout << endl;
cout << "dati valoa rea in care se aproximeaza functia" << endl;
cin >> a;
l=funct(x[1]);
for (k=2;k<=n;k++)
{
s=0;
for (i=1;i<=k;i++)
{
p=1;
for (j=1;j<=k;j++)
if (j !=i)
p=p*(x[i] -x[j]);
s=s+funct(x[i])/p;
}
p1=1;
for (m=1;m<=k -1;m++)
p1=p1*(a -x[m]);
l=l+s*p1;
}
cout << "Val. pol. in " << a << " este " << l << endl;
getche();
}
42
43
4.3. Problema nr. 3
Scrie ți un program în C++ care determină valoarea aproximativă a funcției f în valoarea k
prin metoda polinomului Lagrange, pentru funcția f:R R, f(x)=sin(x)+cos(x).
Gradul polinomului n=3, punctele de tabelare (0, f(0)), (2*PI/3, f(2*PI/3)), (4*PI/3,
f(4*PI/3)), (2*PI, f(2*PI)).
k=5*PI/3=5.236.
După calcule punctele de tabelare sunt: (0, 0), (2.094, 0.366), (4.189, 1.366), (6.283, 1).
Programul calculeaz ă f(k) = f(5.236) = 1.4786
44
4.4. Problema nr. 4
Scrie ți un program în C++ pentru obținerea unui polinom de interpolare pornind de la
următoarele date de intrare: n = 3, x = -1, 0, 1, 2 și y = 1, 2, 0, 3.
Programul sursă C++
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
//#include <alloc.h>
//Aloca dinamic spatiu pentru un vector de dimensiune 'size'
double* alloc_vector(int size)
{
double* vect;
vect = (double*)malloc(sizeof( double)* size);
45
return vect;
};
//Citeste elementele unui vector de dimensiune 'size', alocat anterior cu 'alloc_vector()'
int read_vector(int size, double* vect, char *str)
{
int line;
double t;
for(line = 0 ; line < size ; line++)
{
do{
printf("%s[%d] : ", str, line);
fflush(stdin);
}while(scanf("%lf", &t) != 1);
vect[line] = t;
};
return 0;
};
//Afiseaza elementele unui vector
int print_vector(int size, double* vect, char *str)
{
int line;
printf(" \n");
for(line = 0 ; line < size ; line++)
printf("%s[%d] = %lf \n", str, line, vect[line]);
printf(" \n");
return 0;
};
46
//Elibereaza spatiul alocat unui vector cu 'alloc_vector()'
int free_vector(double* vector)
{
free(vector);
return 0;
};
//Initializeaza elementele unui vector alocat, cu valoare din 'ini_val'
int init_vector(int size, double* vect, double ini_val){
int i;
for(i = 0 ; i < size ; i++)
vect[i] = ini_val;
return 0;
}
void coef_poli(int size, double* x, double* y, double* p)
{
int i, j, k, grad_tem p;
double coef_temp; // pastreaza coeficientul TOTAL de pe linga un polinom partial de la
pasul 'i', adica are valoarea conform formulei din curs 'y[i]/produs(x[i] -x[j]),cu j=0:n si j!=i'
double* temp; // vector temporar, pastreaza coeficientii poli noamelor de la pasul 'i'
if ((temp = alloc_vector(size)) == NULL){
printf(" \nEroare alocare…");
exit(0);
}
double* q; // vector de 'manevre' ce pastreaza coeficientii polinoamelor
'inferioare',obtinute din produse partiale la pasul 'j'
if ((q = alloc_vector(size)) == NULL){
printf(" \nEroare alocare…");
47
free_vector(temp);
exit(0);
}
init_vector(size,p,0.0);
for(i = 0 ; i < size ; i++)
{
grad_temp = 0; // la fiecare pas 'i', gradul polinomului 'temp' normal e 0
init_vector(size,temp,0.0); // adica toti coeficientii lui 'temp' sunt 0.0
temp[0] = 1.0; // in afara de termenul liber, adica normal ca 'temp' nu poate fi identic
nul,pentru ca altfel produsele partiale de polinoame ar da polinoame identic nule
coef_temp = y[i];
if (coef_temp != 0)
{ // daca avem functia nula intr -un punct, atunci polinomul partial la pasul 'i' va
fi identic nul => nu -l luam in con sideratie
for(j = 0 ; j < size ; j++)
if (j == i) continue; // pentru determinarea coeficientilor prin produs,pentru evitarea
impartirii la zero,trebuie ca i!=j,conditie evitate prin 'continue'
else
{
coef_temp /= (x[i] – x[j]);
init_v ector(size,q,0.0); // facem identic nul polinomul de manevre
for(k = 0 ; k <= grad_temp ; k++)
{ // calculul propriu -zis al coeficientilor polinoamelor de produse partiale si
'update -area' coeficientilor polinomului 'temp' de la pasul i
q[k+1] + = temp[k];
if (k == 0) q[k] = temp[k]*( -x[j]);
else q[k] += temp[k]*( -x[j]);
temp[k] = q[k];
}
temp[k] = q[k];
grad_temp++;
}
for(j = 0 ; j < size ; j++)
48
p[j] += coef_temp*temp[j]; // polinomul P fiind definit ca suma de pol inoame de la
pasi 'i', facem exact acest lucru, adica insumarea directa in P
}
}
free_vector(temp);
free_vector(q);
}
int main()
{
int n = 1; // n – gradul polinomului de interpolare SAU indicele maxim al punctelor ce
determina functia
printf(" \nDati n : ");
fflush(stdin);
scanf("%d",&n);
double* p; // vector in care vom pastra coeficientii polinomiali Lagrange
if ((p = alloc_vector(n+1)) == NULL){
printf(" \nEroare alocare…");
exit(0);
}
double* x; // vector pe ntru pastrarea punctelor in care avem cunoscuta functia
if ((x = alloc_vector(n+1)) == NULL){
printf(" \nEroare alocare…");
free_vector(p);
exit(0);
}
double* y; // vector pentru pastrarea valorilor functiei in punctele 'x', in care cunoastem
functia
if ((y = alloc_vector(n+1)) == NULL){
printf(" \nEroare alocare…");
free_vector(p);
free_vector(x);
49
exit(0);
}
printf(" \nDati Xi (punctele in care cunoastem functia) : \n");
read_vector(n+1,x,"X");
printf(" \nDati Yi (valorile functiei in punctele Xi) : \n");
read_vector(n+1,y,"Y");
coef_poli(n+1,x,y,p);
printf(" \nCoeficientii polinomiali Lagrange :");
print_vector(n+1,p,"P");
free_vector(x);
free_vector(y);
free_vector(p);
};
50
Polinomul de interpolare este: 1.333333*x3 – 1.5*x2 – 1.833333*x + 2.
4.5. Problema nr.5
Să se determine diferența divizată de ordinul 3 pentru funcția
1 )(2 x xxf și no durile
de interpolare 1, 3, 4 .
Algoritmul
Algoritmul pentru determinarea diferențelor divizate este
Intrări: k = numărul nodurilor de interpolare
x = vectorul cu nodurile de interpolare
f = funcția din metodă
Ieșiri: dd = diferența divizată
{
dd
0
pentru i
1 : k
{
p
1
pentru j
1 : k
dacă
ij atunci
) ( j ixxp p ;
pxf dd ddi/)( }
}
51
Programul sursă C ++
// Diferenta divizata
#include <stdio.h>
#include < conio.h>
#include <iostream >
using namespace std;
#define max 10
float funct(float t)
{
float f;
f=t*t+t+1;
return f;
}
int main ( )
{
float dd,p;
float x[max+1];
int k,i,j;
cout << "dati ordinul diferentei divizate" << endl; cin >> k;
cout << "dati nodurile" << endl;
for (i=1;i<=k;i++)
{
cout << "x[" << i << "]="; cin >> x[i];
}
cout << endl;
dd=0;
for (i=1;i<=k;i++)
{
p=1;
for (j=1;j<=k;j++)
if (j !=i)
p=p*(x[i] -x[j]);
dd=dd+funct(x[i])/p;
}
52
cout << "Diferente divizata de ordinul " << k << " este " << dd << endl;
getche();
}
53
4.6. Problema nr.6
Să se determine diferența divizată de ordinul 3 pentru funcția
1 )(3 xxxf și no durile
de interpolare 1, 3, 4 în primul nod.
Algoritmul
Algoritmul pentru determinarea diferențelor finite este
Intrări: n = numărul nodurilor de interpolare
k = ordinul diferenței finite
x = vectorul cu nodurile de interpolare
h = pasul
Ieșiri: df = diferența finită
{
dacă k este par atunci df
)(1xf
altfel df
)(1xf
pentru j
1 : k
{
c
1
pentru
1 : j
/)1 ( kcc
dacă
) ( jk este par atunci
df
) (1jh xfcdf
altfel
df
) (1jh xfcdf
}
}
Programul sursă C ++
#include <stdio.h>
#include <conio.h>
#include <iostream.h>
#include <math.h>
using namespace std;
#define max 10
float funct(float t)
{
float f;
54
f=exp( -(t*t));
return f;
}
int main ( )
{
float df,c,h;
float x[max+1];
int n,k,i,j,l;
cout << "dati ordinul diferentei finite" << endl;
cin >> k;
cout << "dati numarul de noduri" << endl;
cin >> n;
cout << "dati pasul" << endl;
cin >> h;
cout << "dati primul nod" << endl;
cin >> x[1];
for (i=2;i<=n;i++)
{
x[i]=x[1]+ (i-1)*h;
}
cout << endl;
if (k%2==0)df=funct(x[1]);
else df= -funct(x[1]);
for (j=1;i<=k;j++)
{
c=1;
for (l=1;l<=j;l++)
c=c*(k -l+1)/l;
if ((k -j)%2==0)df=df+c*funct(x[1]+j*h);
else df=df -c*funct(x[1]+j*h);
}
cout << "Diferenta finita de ordinul " << k << " este " << df << endl;
getche();
}
55
56
4.7. Problema nr. 7
Să se scrie polinomul Newton (varianta progresivă și regresivă) pentru funcția
xexf)(
și nodurile 1, 3, 5, 7.
Rezolvare.
Polinomul lui Newton progresiv este:
)2 )(1()21( )1(!31)1()21( )1(!21)21( )1()1 (!3)1()1 (!2)1()1 (!1)1()1 (!)1()1( )21(
3
0332
0221
0113
13 2
12 1
113
1 14
t ttj fCttj fC tj fC etftftfetkfft L
jjjjjj
jjjkk k
de unde se obține
)2 )(1() 3 3 (61)1() 2(21) ( )21(
7 5 35 3 3
t tt e e eette ee tee et Ln
în timp ce polinomul lui Newton regresiv este:
)2 )(1()21( )1(!31)1()23( )1(!21)25( )1()1 (!3)1()1 (!2)3()1 (!1)5()1 (!) ()7( )27(
3
0332
0221
011 73
13 2
12 1
1173
1 14
t ttj fCttj fC tj fC etftftfetkxff t L
jjjjjj
jjjkkkk
n
de unde
).2 )(1() 3 3 (61)1() 2 (21) ( )27(
7 5 37 5 3 7 5 7
t tte e eette e e te e et Ln
Algoritmul
Algoritmul pentru determinarea valorii polinomului Newton pe noduri echidistante în
„primul” nod este
Intrări: n = numărul nodurilor de interpolare
x = vectorul cu nodurile de interpolare
h = pasul
a = numărul în care se calculează valoarea pol inomului
Ieșiri:
1 = valoarea polinomului în a
{
1
)(1xf
57
pentru k
1 : n – 1
{
p
1
pentru
1 : k
)1 ( ap p
1fc
pentru
1 : k
fc fc
dacă k este par atunci
1s
)(1xf
altfel
1s
)(1xf
pentru j
1 : k
{
1c
pentru
j s :1
s skcc /)1 (
dacă
) ( jk este par atunci
) ( 1 11jh xcfs s
altfel
) ( 1 11jh xcfs s
}
fcps /11 1
}
}
Algoritmul
Algoritmul pentru determinarea valorii polinomului Newton pe noduri echidistante în
„ultimul” nod este
Intrări: n = numărul nodurilor de interpolare
x = vectorul cu nodurile de interpolare
h = pasul
a = numărul în care se calculează valoarea po linomului
Ieșiri:
n = valoarea polinomului în a
{
n
)(nxf
pentru k
1 : n – 1
{
58
p
1
pentru
1 : k
)1 ( ap p
1fc
pentru
1 : k
fc fc
dacă k este par atunci
1s
) ( knxf
altfel
1s
) (knxf
pentru j
1 : k
{
1c
pentru
j s :1
s skcc /)1 (
dacă
) ( jk este par atunci
) ( 1 1 jh xcfs skn
altfel
) ( 1 1 jh xcfs skn
}
fcpsn n /1
}
}
Programul sursă C ++
// Polinomul de interpolare Newton -progresiv
#include <stdio.h>
#include < conio.h>
#include <iostream.h>
#include <math.h>
using namespace std;
#define max 10
float funct(float t)
{
float f;
f=exp( -(t*t));
return f;
59
}
int main ( )
{
float np,df,c,h,p,fc,a;
float x[max+1];
int n,k,i,j,l;
cout << "dati valoarea unde se aproximeaza functia" << endl;
cin >> a;
cout << "dati numarul de noduri" << endl;
cin >> n;
cout << "dati pasul" << endl;
cin >> h;
cout << "dati primul nod" << endl;
cin >> x[1];
for (i=2;i<=n;i++)
{
x[i]=x[1]+(i -1)*h;
}
cout << endl;
np=funct(x[1]);
for (k=1;k<=n -1;k++)
{
fc=1;
for (l=1;l<=k;l++)
fc=fc*l;
p=1;
for (l=1;l<=k;l++)
p=p*(a -l+1);
if (k%2==0)df=funct(x[1]);
else df= -funct(x[1]);
for (j=1;i< =k;j++)
{
c=1;
for (l=1;l<=k;l++)
c=c*(k -l+1)/l;
60
if ((k -j)%2==0)df=df+c*funct(x[1]+j*h);
else df=df -c*funct(x[1]+j*h);
}
np=np+p*df/fc;
}
cout << "Valoarea polinomului in " << a << " este " << np << endl;
getche ();
}
61
Programul sursă C ++
// Polinomul de interpolare Newton -regresiv
#include <stdio.h>
#include <conio.h>
#include <iostream.h>
#include <math.h>
using namespace std;
#define max 10
float funct(float t)
{
float f;
f=exp( -(t*t));
return f;
}
int main ( )
{
float np,df,c,h,p,fc,a;
float x[max+1];
int n,k,i,j,l;
cout << "dati valoarea unde se aproximeaza functia" << endl;
cin >> a;
cout << "dati numarul de noduri" << endl;
cin >> n;
cout << "dati pasul" << endl;
cin >> h;
cout << "dati primul nod" << endl;
cin >> x[1];
for (i=2;i<=n;i++)
{
x[i]=x[1]+(i -1)*h;
}
cout << endl;
np=funct(x[n]);
for (k=1;k<=n -1;k++)
{
62
fc=1;
for (l=1;l<=k;l++)
fc=fc*l;
p=1;
for (l=1;l<=k;l++)
p=p*(a+l -1);
if (k%2==0)df=funct(x[n -k]);
else df= -funct(x[n -k]);
for (j=1;i<=k;j++)
{
c=1;
for (l=1;l<=k;l++)
c=c*(k -l+1)/l;
if ((k -j)%2==0)df=df+c*funct(x[n -k]+j*h);
else df=df -c*funct(x[n -k]+j*h);
}
np=np+ p*df/fc;
}
cout << "Valoarea polinomului in " << a << " este " << np << endl;
getche();
}
63
64
5. BIBLIOGRAFIE
[1]. Bucur, C.M., Metode numerice , Ed. Facla, Timișoara, 1973.
[2]. Coman, G., Analiză numerică , Ed. Libris, Cluj, 1995.
[3]. Dodescu, Gh., Toma, M., Metode de calcul numeric , E. D. P., București, 1976.
[4]. Dodescu, Gh., Metode numerice în algebră , Ed. tehnică, București, 1979.
[5]. Ichim, I., Marinescu, G., Metode de aproximare numerică , Ed. Academiei R. S. R.,
București, 1986.
[6]. Ignat, C., Ilioi, C., Jucan, T., Elemente de informatică și calcul numeric , Univ . „Al.
I. Cuza”, Iași, Fac . de Matematică, 1989.
[7]. Scheiber, E., Metode numerice , Univ . Transilvania din B rașov, Fac ultatea de
Matem atică – Informatică, (electronic).
[8]. Toma, M., Odăgescu, I., Metode numerice și subrutine , Ed. Tehnică, București, 1980.
[9]. Talmaciu , M., Patriciu , A.M., Calcul numeric , Editura Pim, Iași, 2002.
[10]. Talmaciu , M., Patriciu , A.M., Analiză numerică cu exemple în limbajul C și
Pascal, Editura EduSoft, Bacău, 2005.
[11]. Talmaciu ,M., Patriciu , A.M., Recipiente numerice în Pascal și C, Editura
EduSoft, Bacău, 2005.
[12]. Talmaciu , M., Patriciu , A.M., Metode numerice, Volumul I, Editura Tehnica
Info, Chișinău, 2005.
Copyright Notice
© Licențiada.org respectă drepturile de proprietate intelectuală și așteaptă ca toți utilizatorii să facă același lucru. Dacă consideri că un conținut de pe site încalcă drepturile tale de autor, te rugăm să trimiți o notificare DMCA.
Acest articol: Str. Calea Mărășești, nr. 157, Bacău, 600115 [608895] (ID: 608895)
Dacă considerați că acest conținut vă încalcă drepturile de autor, vă rugăm să depuneți o cerere pe pagina noastră Copyright Takedown.
