Transcript MNMSD03

CURS 3 APROXIMAREA FUNCTIILOR DE O VARIABILA (II)
b) Interpolarea cu functii spline
Functiile spline constituie un instrument matematic de data relativ recenta, dar cu o larga
aplicabilitate in procesarea pe calculator a datelor numerice. Termenul “spline”, de origine
engleza, a fost adoptat ca atare in literatura de specialitate si reprezinta denumirea unei rigle
flexibile prevazuta cu o serie de greutati de-a lungul ei, folosita de catre desenatori in trasarea
unor curbe care trec printr-o serie de puncte.
Functia poligonala este unul din cele mai simple exemple de functii spline.
Functiile spline sunt caracterizate prin formele lor pe
subintervale dintre douã noduri (care pot fi diverse functii
cunoscute) si prin anumite conditii de racordare in noduri.
f(b)
s(x)
f(a)
f(x)
x1 = a
xn = b
A doua conditie contine conditia
de racordare in noduri:
(k)
(k)
sm,i(xi+1) = sm,i+1(xi+1) k = 0, 1,…, m-1
Fie [a,b]  R un interval din dreapta reala si xi o retea
de puncte de diviziune (x1 = a, xn = b). Se noteaza
cu Ii subintervalele [xi, xi+1).
Prin definitie, functia s; [a,b]  R se numeste functie
spline polinomiala de ordinul m dacã:
a) Restrictiile acesteia pe subintervalele Ii sunt
polinoame de gradul m;
b) s este derivabilã de (m-1) ori pe intervalul (a,b).
b1) Functia spline de ordinul 1 (linia poligonala)
Este functia cea mai simpla, formate din segmente de dreapta de ecuatie:
yi+1 - yi
s1,i(x) = mix +ni , x  [xi, xi+1), cu panta pe subintervalul Ii: mi =
.
hi
In relatia respectiva hi = xi+1 – xi este marimea subdiviziunii. Cum s1,i(x) = yi se obtine ni = yi – mixi.
Asadar, s1,i(x) = yi + mi(x – xi) , x  [xi, xi+1), i = 1,…, n-1
b2) Functia spline de ordinul 2
Este formata din segmente de parabola racordate in noduri pana la derivata de ordinul intai
inclusiv, adica pentru doua intervale alaturate, atat valorile celor doua functii spline cat si
valorile derivatelor lor de ordinul intai, in nodul comun, sunt egale.
s2,i(x) = yi + mi(x – xi) + ai(x – xi)2 , x  [xi, xi+1), i = 1,…, n-1
‘ (xi+1), rezultã urmatoarele valori pentru coeficientii ai:
‘ (xi+1) = s2,i+1
Cum s2,i(xi+1) = yi+1 si s2,i
ai =
yi+1 - yi
2
hi
m
- i
hi
yi+1 - yi
mi + mi+1 = 2
hi
pentru pantele in noduri.
i = 1,…, n-1
b3) Functia spline de ordinul 3 (functia spline cubica)
Este una din cele mai utilizate functii ce are derivate continue pana la ordinul doi inclusiv,
permitand astfel si calculul razei de curbura.
Aceasta are expresia: s3,i(x) = yi + mi(x – xi) + bi(x – xi)2 + ai(x – xi)3 , x  [xi, xi+1), i = 1,…, n-1
Conditiile de continuitate in noduri sunt:
‘ (x ) si s “ (x ) = s “ (x ) rezultand urmatoarele valori pentru
‘ (xi+1) = s3,i+1
s3,i(xi+1) = yi+1 , s3,i
i+1
3,i i+1
3,i+1 i+1
coeficientii ai si bi:
ai =
mi+1 + mi
bi = 3
h2i
yi+1 - yi
hi2
-2
+
yi+1 - yi
hi3
mi+1 + 2mi
i = 1,…, n-1
3
hi
Din a treia conditie se obtine sistemul care dã pantele pe noduri mi:
imi-1 + 2mi + imi+1 = di , in care:
i =
di =
hi
hi-1 + hi
i =
hi-1
hi-1 + hi
y -y
y –y
3
(hi-1 i+1 i + hi i i-1 )
hi-1 + hi
hi
hi-1
Intrucat relatia in mi-1, mi, mi+1 reprezinta (n-2) ecuatii cu n necunoscute, sistemul trebuie
completat cu incã doua conditii pentru a suplini conditiile de racordare in prima si a doua derivatã
ce nu mai pot fi scrise pentru nodul xn (spre exemplu se pot preciza pantele m1 si mn).
APLICATIA 1
(EXEMPLU DE CALCUL ANALITIC)
Se considerã urmãtoarele citiri:
i
1
2
3
xi
0
2
3
yi
1
2
0
Sã se aproximeze, utilizand interpolarea cu functii spline cubice, valoarea lui f(2.5).
REZOLVARE
Se aplica relatiile: i =
di =
hi
hi-1 + hi
i =
hi-1
hi-1 + hi
pentru i = 2 , obtinandu-se:
y -y
y –y
3
(hi-1 i+1 i + hi i i-1 )
hi-1 + hi
hi
hi-1
45o
1
2
; 2 =
; d2 = - 7
3
3
2
Se impun pantele la capete m1 = 1 si m3 = 0, iar din relatia: imi-1 + 2mi + imi+1 = di
rezultã: m2 = - 23
12
y3 – y2
m3 + m 2
y2 – y1
m2 + m 1
2
a
=
-2
a1 =
2
3 = 2.0833
2
3 = - 0.4792
2
h
h
h2
2
h1
1
y3 – y2 m3 + 2m2
y2 – y1 m2 + 2m1
+
= - 2.1667
b2 = 3
+
= 0.7083
b1 = 3 2
2
3
3
h2
h2
h1
h1
h1 = x2 – x1 = 2 ; h2 = x3 – x2 = 1 ; 2 =
Inlocuind aceste valori in relatia: s3,i(x) = yi + mi(x – xi) + bi(x – xi)2 + ai(x – xi)3 , rezulta
expresiile functiei spline pe subintervalele [0,2] si [2,3]:
s3,1(x) = y1 + m1(x – x1) + b1(x – x1)2 + a1(x – x1)3 = 1 + x + 0.7083x2 – 0.4792x3
s3,2(x) = y2 + m2(x – x2) + b2(x – x2)2 + a2(x – x2)3 = - 19.5 + 31.75x – 14.66x2 + 2.083x3
Cum x = 2.5  [2,3], s3,2(2.5) = 0.7608
c) Aproximarea prin metoda celor mai mici patrate
Se stie ca datele obtinute prin masurarea cantitativa a variabilelor din fenomenele fizice
studiate sunt afectate de o serie de erori: erori de observatie, erori de masurare sau erori de
inregistrare. In aceste situatii aproximarea functiei printr-un polinom de interpolare nu este
recomandata, deoarece polinomul se suprapune exact in toate punctele determinate
experimental si prin urmare erorile continute in datele masurate se vor pastra prin
interpolare. Aproximarea functiilor prin metoda celor mai mici patrate conduce la
determinarea unei relatii functionale care elimina intr-o oarecare masura erorile. Criteriul de
aproximare il reprezinta minimizarea sumei:
m
S=
[y – g(x )]

i=1
i
i
2
= minim
n
valorile xi si yi fiind date. Aproximanta g(x) se scrie sub forma: g(x) =  akgk(x)
k=1
gk(x) sunt functii cunoscute, liniar independente, iar ak – parametri nedeterminati.
k = 1, 2…, n
c1) Aproximarea liniara a functiilor prin metoda celor mai mici patrate
In acest caz polinomul g(x) se ia de forma:
g(x) = P1(x) = a1 + a2x
Problema revine la minimizarea functiei:
m
S(a1, a2) =
[y – (a +a x )]

i=1
i
1
2 i
2
in raport cu parametrii a1 si a2.
Pentru aceasta este necesar ca derivatele partiale sa se anuleze:
S
si
a1 = 0
S
a2 = 0
Acest sistem este echivalent cu:
m
m
a1m + a2 xi = yi
m
i=1
m
i=1
m
a1x + a2 x2 = x y
i
i i
i
i=1
i=1
i=1
si permite determinarea unica a coeficientilor a1 si a2 daca se cunosc valorile (xi, yi), i = 1,…, m,
unde m reprezinta numarul punctelor citite (numarul masuratorilor).
c2) Aproximarea functiilor printr-un polinom de grad (n-1) utilizand metoda celor mai
mici patrate
Judecand intr-un mod similar aproximarii liniare expuse, se obtine cazul aproximarii functiei
necunoscute printr-un polinom de gradul (n-1), cu n < m:
n
g(x) = Pn-1(x) =  akxk-1 = a1 + a2x + a3x2 + …+ anxn-1
k =1
sistem cu n ecuatii si n necunoscute, exprimat matriceal prin urmatoarea ecuatie:
m
xi

i =1
m
m
m
i =1
m
…

i =1
…
…
…
m
m

i =1
…
m
n-1

2
xi
…
2
xi
xi

i =1
xi

i =1
m

i =1
n
xi
3
xi

i =1
n+1
xi
…
permitandu-se determinarea coeficientilor ak.
m
xi

i =1
m
n-1

i =1
n
xi
…
m

i =1
2n-2
xi
m
yi

i =1
a1
a2
……
an
m
=
xiyi

i =1
…
m n-1
xi yi

i =1
(EXEMPLU DE CALCUL ANALITIC)
APLICATIA 2
Se considera tabelul cu puncte date, generate de o functie necunoscuta:
i
1
2
3
4
5
xi
0
0.25
0.5
0.75
1
yi
1
1.2840
1.6847
2.1170
2.7183
Sã se aproximeze aceste date cu un polinom de gradul doi, utilizand metoda celor mai
mici patrate.
REZOLVARE
In cazul de fatã avand 5 citiri  m = 5, iar n = 3. Polinomul de gradul doi este de forma:
P2(x) = a1 + a2x + a3x2. Se aplica sistemul liniar prezentat anterior, care se particularizeaza
pentru n = 3. Se va obtine asadar sistemul:
5
5
i =1
i =1
5a1 + a2  xi + a3 
5
5
i =1
5
i =1
5
 a1  xi + a2 
a1 
2
xi
i =1
+ a 2
2
xi
3
xi
i =1
2
xi
5
5
= yi
5a1 + 2.5a2 + 1.875a3 = 8.768
care se mai scrie:  2.5a1 + 1.875a2 + 1.5625a3 = 5.4514
1.875a1 + 1.5625a2 + 1.3828a3 = 4.4015
i =1
+ a3 
3
xi
+ a3 
4
xi
i =1
5
i =1
5
=  xiyi
i =1
5 2
=  xi yi
i =1
Rezulta: a1 = 1.0052, a2 = 0.8641, a3 = 0.8437
Polinomul va fi: P2(x) = 0.8437x2 + 0.8641x + 1.0052
5
Eroarea va fi:
i=1[y – P (x )]
i
2
i
2
= 0.0003
PROCEDURA MATLAB PENTRU GASIREA COEFICIENTILOR UNUI
POLINOM UTILIZAND METODA CELOR MAI MICI PATRATE
Funcţia polyfit găseşte coeficienţii unui polinom (o curbă) care aproximează un set de
date în sensul algoritmului celor mai mici pătrate:
p = polyfit(x,y,n)
x şi y sunt vectorii care conţin setul de date iar n este ordinul polinomului ai cărui
coeficienţi vor fi furnizaţi la apelarea funcţiei.
COMANDA
>> x = [0, 0.25, 0.5, 0.75, 1]
>> x = [0, 0.25, 0.5, 0.75, 1]
>> y = [1, 1.2840, 1.6487, 2.1170, 2.7183]
>> y = [1, 1.2840, 1.6487, 2.1170, 2.7183]
>> p = polyfit (x,y,2)
>> p = polyfit (x,y,3)
(polinom de gradul trei)
(polinom de gradul doi)
REZULTATE OBTINUTE
p=
0.8437
p=
0.8641
1.0052
0.2789
0.4253
1.0141
0.9999
2.8
2.6
2.4
2.2
2
1.8
1.6
1.4
1.2
1
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Comanda
suplimentara
care se da daca se
doreste ca punctele prin
care trece polinomul de
interpolare sa apara pe
desen
>> t = 0: .25 :1;
Comanda in MATLAB
pentru trasarea graficului
polinomului de interpolare
de gradul doi
>> y = 0.8437*t.^2 + 0.8641*t+1.0052;
>> plot (t,y)
>> grid on
>> plot (t, y, ’-r’ ,t ,y ,’ok’)