Aproximarea functiilor prin interpolare: Polinomul de interpolare Lagrange, Newton

Aproximarea functiilor prin interpolare

Tematica lucrarii: 1. Generalitati;

2. Polinomul de interpolare Lagrange;

3. Polinomul de interpolare Newton;

4. Interpolarea cu functii Spline;

5. Concluzii si intrebari.

1. Generalitati

Fie o functie reala de variabila reala, data print-o multime de valori discrete asociate unei multimi de valori ale argumentului. Prin procesul de interpolare pot fi determinate valori ale functiei corespunzatoare unor valori oarecare ale argumentului care se gasesc intre valorile date. Astfel o functie discreta, sau data prin valori discrete, poate fi aproximata printr-o functie de interpolare de o anumita forma.

Functia de interpolare poate fi definita pentru tot domeniul corespunzator multimii de valori discrete date, sau poate fi definita pe portiuni. Functiile de interpolare pot avea diferite forme, apartinand unor anumite clase de functii. Acestea ofera diferite avantaje si dezavantaje din punctul de vedere al erorilor de aproximare locala si globala. De asemenea functia de interpolare poate avea diferite proprietati din punctul de vedere al gradului de netezime, adica rangul pana la care contine derivate continue.

2. Polinomul de interpolare Lagrange

Interpolarea polinomiala utilizeaza ca functie de interpolare un polinom de un anumit grad. Gradul acestuia depinde de numarul punctelor date. Odata cu cresterea gradul polinomului de interpolare creste si eroarea de aproximare in punctele departate de acestea.

Coeficientii polinomului Lagrange sunt determinati pe baza conditiilor de coincidenta intre valorile functiei de interpolat si polinomul de interpolare in punctele multimii discrete date. Unul dintre avantajele interpolarii Lgrange este posibilitatea utilizarii unei multimi de puncte neechidistante.

Fie o diviziune cu n+1 puncte a intervalului [a,b]: x , x , x , , x_n, cu x = a si x_n = b in care sunt date valorile functiei reale definite pe intervalul [a,b]:

y = f(x y = f(x ) ; y = f(x ) ; ; y_n = f(x_n) ; y_n = f(x_n (1)

Polinoamele Lagrange, notate L_n(x), au forma:

(2)

si sunt definite de conditiile:

L_n(x_i) = y_i pentru i = 0,1,2,,n (3)

Polinoamele elementare p_i(x) sunt definite cu ajutorul simbolului lui Kroneker _ij astfel:

p_i(x) = _ij daca x = x_i si (4)

= 0 daca x = x_j , j i

Se poate observa ca polinomul p_i(x) are ca radacini punctele x ,x ,,x_i-1,x_i+1,,x_n , iar in punctul x_i trebuie sa ia valoarea 1. Daca acestea sunt scrise sub forma:

p_i(x) = c_i(x-x )(x-x ) (x-x_i-1)(x-x_i+1) (x-x_n (5)

si determinam coeficientii c_i din conditiile: p_i(x_i) = 1 pentru diferite valori ale lui i, atunci rezulta:

(6)

si respectiv:

(7)

Daca utilizam notatiile:

ⁱ_n+ (x) = (x-x )(x-x ) (x-x_i-1)(x-x_i+1) (x-x_n (8)

atunci polinomul Lagrange poate fi scris sub forma condensata:

(9)

Din relatia de definitie (5.32) se poate observa ca in punctele x_i polinoamele elementare p_i(x_j) au valoarea unu, iar polinomul de interpolare Lagrange are valori egale cu valorile functiei y_i . Deci prin aproximarea cu polinoame Lagrange esantioanele functiei sunt conservate, iar in intervalele dintre puntele date avem o variatie de tip polinomial. De aceea in cazul polinoamelor de grad superior ne putem astepta la o eroare de aproximatie importanta in zona de mijloc a intervalelor de interpolare.

Programarea in Mathcad

Fie o multime discreta care contine N puncte. De exemplu N := 9.

Alegem doua variabile sir, utilizate ca indici:

i := 0 .. N si j := 0 .. N (10)

Notam cu variabila indexata x_i vectorul punctelor diviziunii si cu y_i vectorul valorilor functiei date. Sa consideram ca x_i parcurge valorile lui i, iar functia data are o evolutie oarecare, generata aleator de functia predefinita rnd(1) unde valoarea 1 indica capatul domeniului de valori generate incepand cu valoarea zero:

y_i := rnd(1) (11)

Sa alegem pasul s := 1 pentru variatia valorilor din domeniul de definitie a polinomului de interpolare Lagrange. Consideram pentru capatul acestui interval valoarea N. Alegem variabila indexata z_k pentru vectorul valorilor argumentului polinomului de interpolare Lagrange. Rezulta atunci:

k := 0.. N/s si z_k := s k (12)

Tinand cont de teoria polinomului de interpolare Lagrange rezulta ca acesta poate fi exprimat in limbaj Mathcad prin ecuatia recursiva:

(

Problema de calcul

Se cere implementarea si rezolvarea in limbaj Mathcad a problemei de interpolare, utilizand polinomul Lagrange, pentru functia y(x). In acest scop se vor trasa curbele de variatie pentru functia data si pentru polinomul de interpolare pe acelasi grafic. Sa se urmareasca apoi eroarea de aproximare pe tot domeniul de variatie pentru cele doua functii.

Fig.1 Variatia functiei discrete y(x) (interpolata liniar).

Fig.2 Variatia polinomului Lagrange si a functiei y(x)

Se cere:

Implementarea in limbaj MathCad a algoritmului de interpolare cu polinoame Lagrange;

Studiul modului de alegere a punctelor si a numarului de puncte, pentru o diviziune a intervalului de variatie, pentru o anumita functie data;

Studiul calitatii aproximarii pentru diferite tipuri de functii.

3. Polinomul de interpolare Newton

Formulele de interpolare ale lui Newton se aplica in cazul distributiilor de valori ale functiei corespunzatoare unor puncte echidistante (diviziune uniforma). Pentru exprimarea concisa a formulelor de interpolare cu puncte echidistante vom itroduce notiunea de diferenta finita, si putere generalizata.

Daca notam cu y valoarea functie f in punctul x: y = f(x) si cu x notam cresterea argumentului x pe parcursul unui pas de discretizare h, atunci diferenta finita de ordinul unu a functiei f in punctul x, notata cu y, este definita de relatia:

y = f(x) = f(x+x) - f(x)      (14)

Analog notam cu y diferenta finita de ordinul doi, a functiei f in punctul x :

y = [f(x+x) - f(x)] = f(x+x) - f(x) =

[f(x+2x) - f(x+x)] - [f(x+x) - f(x)] = f(x+2x) - f(x)      (15)

Diferenta finita de ordin n o definim prin relatia recursiva:

ⁿy ^n-1y) (16)

Notiunea de putere generalizata definita de relatia:

x^[n] = x[x-h][x-2h] [x-(n-1)h]      (17)

unde h este o constanta fixa, in cazul nostru reprezentand pasul diviziunii uniforme. Prin definitie x^[ = 1, iar daca h = 0 puterea generalizata este identica cu puterea propriu-zisa.

Polinomul de interpolare Newton de ordinul n de speta I are forma:

P_n(x) = a + a (x-x ) + a (x-x x-x ) + + a_n(x-x )(x-x ) (x-x_n-1 (18)

Daca utilizam notatia puterii generalizate polinomul se poate scrie:

P_n(x) = a + a (x-x + a (x-x + + a_n(x-x ^[n]     

Pentru determinarea coeficientilor a , a , ,a_n punem conditiile:

P_n(x ) = y (20)

P_n(x y

P_n(x y

ⁿP_n x ⁿy

Rezolvand sistemului de ecuatii rezulta coeficientii a , a , ,a_n si polinomul de interpolare Newton de speta I devine:

(21)

Daca in aceasta relatie facem substitutia:

(22)

rezulta expresia de forma:

(23)

Caracteristici ale polinoamelor de interpolare Newton

Coeficientii polinomului Newton de speta I sunt determinati pe baza coincidentei diferentelor finite de diferite ordine, ale functiei aproximate si respectiv ale polinomului de interpolare, avand ca punct de referinta capatul din stanga al intervalului pentru care sunt date valorile discrete ale functiei. Acest mod de lucru duce la doua proprietati importante, si anume, obtinerea unei foarte bune aproximari a modului de variatie a functiei in vecinatatea capatului de interval, pe de o parte si posibilitatea extrapolarii catre dreapta pentru anumite clase de functii. Astfel panta unei functii si modul ei de variatie este foarte bine aproximata prin metoda lui Newton.

Programarea in MathCad

Pentru a studia modul de implementare in MathCad a interpolarii functiilor prin metoda lui Newton consideram exemplul de mai jos. Ne propunem sa aproximam functia exponentiala:

f(x) = e^x (24)

Consideram ca functia este data prin valori discrete pe un interval restrans:

(a,b) = (3 ; 3,7) (25)

Ne propunem ca prin metoda lui Newton sa realizam iterpolarea functiei pe intervalul dat si apoi sa realizam extrapolarea acesteia pe o extensie a intervalului:

(b,c) = ( 5) (26)

Daca consideram o diviziune cu patru segmente a intervalului dat (N = 4) rezulta pasul si nodurile diviziunii:

k := 0..N x_k := a + hk (27)

Luam pentru valorile discrete ale functiei de aproximat valorile corespunzatoare nodurilor diviziunii:

fd_k := f(x_k (28)

Daca consideram ca vectorul diferentelor finite de ordinul zero este definit chiar de vectorul valorilor functiei in punctele diviziunii:

(29)

Atunci pentru calculul diferentelor finite pentru ordinele 1..4 ale functiei f(x) putem utiliza urmatoarea ecuatie vectoriala recursiva: (30)

Reamintim ca exponentul <k> plasat intre paranteze ascutite indica vectorul corespunzator pasului de calcul k, iar indicele inferior indica numarul componentei vectorului (numaratoarea facandu-se incepand cu zero. Ecuatia recursiva de mai sus determina diferentele finite, din aproape in aproape, obtinand valori succesive pe baza celor deja calculate. Fiecare diferenta fiunita va fi calculata prin diferenta diferentelor finite precedente, avand ca referinta totdeauna capatul din stanga al subintervalului. Astfel la fiecare pas (ordin al diferentei finite) numarul de valori nenule se micsoreaza cu cate o unitate.

Pentru determinarea coeficientilor polinomului lui Newton pot fi utilizate doua forme de reprezentare. Prima se bazeaza pe substitutia cu variabila q , iar cealalta pe notiunea de putere generalizata. In limbaj MathCad noua variabila q va fi definita ca functie pentru a pastra dependenta de x. Tinand cont ca a = x este capatul din stanga al intervalului, rezulta:

(31)

In MathCad polinomul Newton de speta I va fi definit ca o functie determinata de urmatoarea relatie:

(32)

unde fd reprezinta valoarea functiei in capatul de referinta al intervalului;

[df^<k> reprezinta prima componenta a vectorului diferentelor finite de ordinul k;

k si m sunt variabile sir intregi, definite inaintea expresiei de mai sus:

k := 1..3 ; m := 1..3 (33)

Expresia conditionala care apare in cadrul produsului indexat dupa m este utilizata pentru a putea genera produse cu un numar diferit de factori pentru fiecare pas k.

A doua modalitate de a construi polinomul Newton de speta I este folosirea notiunii de putere generalizata. La implementarea in cadrul MathCad, pentru usurinta exprimarii vom modifica expresia puterii generalizate astfel incat aceasta sa includa numitorul dat de produsul intre factorial si puterea pasului h. Astfel expresia modificata a puterii generalizate de ordinul k va fi data de expresia:

(34)

Rezulta polinomul lui Newton de speta I exprimat in functie de puterea generalizata a lui

(x-a):

(35)

Forma numerica a implementarii polinomului Newton poate fi usor reprezentata grafic prin facilitatile specifice interpretorului MathCad. Sa alegem intervalul (3 ; 5) pentru interpolarea, respectiv extrapolarea functiei f(x) = e^x si sa construim urmatorul sir de valori pentru variabila independenta x :

a b := 5 N := 50 k := 0..N (36)

Sirul de 50 de valori are o rezolutie suficienta pentru a putea trasa si verifica grafic modul de aproximare a functiei continue f(x) = e^x prin interpolare si respectiv extrapolare, utilizand cateva valori discrete si polinomul Newton de speta I.

Fig.3. Variatia polinomului Newton de sp.I si functiei e^x.

Se cere:

Implementarea in limbaj MathCad a algoritmului de interpolare cu polinoame Newton;

Studiul modului de alegere a punctelor si a numarului de puncte, pentru o      diviziune a intervalului de variatie, pentru o anumita functie data;

Studiul calitatii aproximarii pentru diferite tipuri de functii.

4. Interpolarea cu functii Spline

Interpolarea cu functii Spline consta de regula intr-o interpolare polinomiala pe portiuni, realizand insa cuplajul neted intre segmentele de curba cu gradul de continuitate cerut (pentru functie si derivatele acesteia). Interpolarea cu functii Spline, pe de o parte elimina dezavantajul mentionat in cazul interpolarii polinomiale (erori mari de aproximare daca numarul punctelor date este mare), iar pe de alta parte respecta cerintele de continuitate cerute pentru functie si derivatele ei. Interpolarea pe subintervale poate fi realizata si cu alte tipuri de functii ca cele trigonometrice sau cele din clasa functiilor speciale. Printre cele mai utilizate sunt functiile Spline cubice care utilizeaza pentru interpolare polinoame de gradul trei.

Fie f(x) o functie reala de variabila reala definita pe intervalul [a,b] : f : [a,b] --> R , o diviziune a intervalului: a = x , x , x , , x_n = b , pentru care cunoastem valorile esantioanelor in punctele diviziunii: f(x_i) = y_i , i = 1, 2, , n .

Proprietatile functiilor Spline cubice:

a) multimea functiilor Spline notata cu S apartine clasei de functii cu derivate de ordinul doi continui pe intreg intervalul [a,b]:

S C²[a,b] (37)

b) in punctele diviziunii functiile S(x) sunt egale cu esantioanele functiei interpolate f(x):

f(x_i) = S(x_i pentru i = 1, 2, , n (38)

c) derivatele de ordinul doi ale functiilor Spline au norma 'patrat integral' minima :

pentru orice functie g(x) (39)

Realizarea efectiva a interpolarii

Functiile Spline cubice utilizeaza polinoame de gradul trei pe portiuni. Rezulta ca derivatele de ordinul doi au variatie liniara pe portiuni.

Pentru interpolarea Spline se impun cerintele urmatoare.

* Derivata a doua a functiei globale de interpolare pe tot intervalul sa fie continua

(S C²[a,b]). Deci grafic aceasteia va consta intr-o linie franta;

* Limitele laterale ale derivatei de ordinul intai sa fie egale in punctele diviziunii intervalului de interpolare (conditie de netezime):

lim_{x -> xi ; x>xi} S'(x) = lim_{x -> xi ; x<xi} S'(x) (40)

* in punctele diviziunii functiile Spline S(x) sunt egale cu esantioanele functiei interpolate f(x) (conditie de coincidenta):

f(x_i) = S(x_i pentru i = 1, 2, , n (41)

Realizarea efectiva a interpolarii cu functii Spline se conformeaza etapelor urmatoare:

Determinarea valorilor derivatelor de ordinul doi in punctele diviziunii intervalului prin rezolvarea unui sistem de ecuatii construit pe baza conditiei de netezime (37).

Determinarea functiilor Spline, pe fiecare portiune, utilizand valorile derivatelor de ordinul doi pe baza conditiei de coincidenta (41).

Programarea in Mathcad

In cadrul programului Mathcad avem la dispozitie urmatoarele functii predefinite:

cspline(x,y) pentru calculul vectorului corespunzator derivatelor de ordinul doi, avand ca argumente vectorul x al punctelor diviziunii si vectorul y al valorilor discrete ale functiei de interpolat;

interp(w,x,y,z) pentru calculul functiei de interpolare Spline in punctele definite de vectorul z. Prin vectorul w vor fi transmise valorile corespunzatoare functiei cspline.

Pentru calculul valorilor functiei de interpolare Spline in punctele corespunzatoare coordonatelor vectorului z utilizam comenzile Mathcad:

w := cspline(x,y) (42)

S_k := interp(w,x,y,z_k (43)

Problema de calcul Calculati si trasati grafic in Mathcad variatia functiei de interpolare Spline, pentru functia y(x) pe intervalul [0,N].

Trasati pe acelasi grafic variatiile functiei interpolate si functiei de interpolare si urmariti eroarea de aproximare.

Fig.4 Interpolarea cu functii Spline cubice

Se cere:

* Realizarea interpolarii cu functii Spline, utilizand instructiunile MathCad;

* Studiul comparativ a interpolarii cu functii Spline si polinoame Lagrange.

5. Concluzii si intrebari

Din modul de definire si experimentele numerice deduceti avantajele si dezavantajele metodei de interpolare cu polinoame Lagrange si cu functii Spline;

Dati exemple de tipuri de functii pentru care este potrivita interpolarea Lagrange sau Spline;

Prezentati cateva din aplicatiile interpolarii functiilor in domeniul preocuparilor matematice sau stiintifice, in general.