Page 1
Seminar 3. Interpolare
Responsabil :
Mihaela Vasile ([email protected] )
Cosmin-Stefan Stoica ([email protected] )
Obiective
În urma parcurgerii acestui seminar studentul va fi capabil să:
folosească diferite metode de interpolare
cunoască diferențele între aceste metode
prezinte avantajele și dezavantajele acestor metode
Breviar teoretic
O funcţie continuă f:[a,b]R, cunoscută numai în punctele x0,x1,…,xn -suportul
interpolării, prin valorile f(x0),f(x1),…,f(xn) este aproximată printr-un polinom
generalizat de interpolare : Pn(x)=a0u0(x)+a1u1(x)…+anun(x),
unde funcţiile liniar independente u0(x),u1(x),…,un(x) formează baza interpolării.
Coeficienţii a0,a1,…,an se determină impunând condiţiile de interpolare:
Pn(xi)=f(xi), i=0:n.
Condiţiile de interpolare conduc la sistemul de ecuaţii liniare:
n:0i,xfxua i
n
0k
ikk
.(un sistem rău condiționat)
1.Interpolare polinomială
1.1.Polinom de interpolare Lagrange
Pentru baza polinomială uk(x)=xk, k=0:n, polinomul de interpolare se calculează cu
formula Lagrange:
n
kii
ik
n
kii
in
k
kn
xx
xx
xfxP
,0
,0
0
.
1.2. Identitatea lui Newton
Page 2
Numim diferenţe divizate ale unei funcţii f în raport cu suportul x0,…,xn:
F0[xi]=f(xi), i=0:n
k
kkkkkk
xx
xxFxxFxxxF
0
1110110
,,,,,,
, k=1:n.
k
ik
ijj
ji
i
kk
xx
xfxxxF
0
,0
10
)(
)(,,, .
Exprimarea polinomului de interpolare cu diferenţe divizate se numeşte formula lui
Newton:
Pn(x)=f(x0)+(x-x0)F1[x0,x1]+(x-x0)(x-x1)F2[x0,x1,x2]+…
+(x-x0)…(x-xn-1)Fn[x0,x1,…,xn]
Eroarea interpolării are expresia:
!1
|)(||)x-(x)x-(x|]x,,x[x,)Fx-(x)x-(x
)1(
n0n01nn0
n
fxE
n
n
, [x0,xn].
Pentru un suport de interpolare echidistant(xi-xi-1=h, xi=x0+ih) se definesc diferenţe
nedivizate:
- progresive
i1ii1iiiii ffxfxfxfhxffxf
i
1k
1i
1k
i
kfff
- regresive
1ii1iiiiii ffxfxfhxfxffxf
1i
1k
i
1k
i
kfff
- centrate
2
1i
2
1i
2
1i
2
1i
iiii
ffxfxf
2
hxf
2
hxffxf
2
1i
1k
2
1i
1k
i
kfff
Intre diferenţele nedivizate şi cele divizate există relaţiile:
niin
n
2
ni
n
ni
n
i
nx,,xFh!nfff
Page 3
Formulele de interpolare Newton-Gregory:
0
n
0
2
001 fn
uf
2
uf
1
ufup
,
x=x0+uh, 0<u<1
n
nn
n
2
nn2 fn
u1f
2
uf
1
ufup
,
x=xn-uh, 0<u<1
1.3. Polinoame ortogonale
Def: Un polinom este ortogonal <pi,pj> =0; oricare i≠j si ||pi||=1.
b
aijjiji dxxwpppp )(,
Exemple de polinoame ortogonale:
1) Cebasev Tn+1-2xTn+Tn-1=0; T0=1; T1=x;
2) Legendre (n+1)Ln+1-(2n+1)xLn +nLn-1=0; L0=1; L1=x;
3) Laguerre Gn+1-(2n+1-x)Gn+n2Gn-1=0; G0=1; G1=1-x;
4) Hermite Hn+1-2xHn+2nHn-1=0; H0=1; H1=2x;
Proprietăți
1) Orice polinom ortogonal are rădăcinile în [a,b] reale și distincte.
2) Orice polinom ortogonal este ortogonal cu orice polinom neortogonal de grad mai
mic decât el.
Demonstrație P2)
0=<pn,pj>; j<n;
<pn,pk>=0=><pn,
k
i
i
i xa1
> =0 =>
k
i 1
< pn,aixi>=0 => <pn,aix
i>=0 => este ortogonal
1.4. Interpolare cu polinoame Cebâsev
Baza interpolării este:
xT,,xT,xT,2
1n21 .
Polinomul generalizat de interpolare este:
xTaxTa2
axP nn11
0
n
Suportul interpolării este constituit de rădăcinile polinomului Cebâşev de grad n+1:
Tn+1(xk)=0, n:0k,2n2
1k2cosxk
Sistemul generat de condiţiile de interpolare este ortogonal şi are soluţia:
Page 4
n
0k
n
0k
k02n2
1k2cosf
1n
2xf
1n
2a
n
0k
jk
n
0k
kjkj2n2
1k2cosTxf
1n
2xTxf
1n
2a
1.5.Interpolare trigonometrică
Considerăm funcţia f periodică cu perioada 2. Suportul interpolării este reprezentat de
2n+1 puncte echidistante în [0, 2]:
n2:0k,1n2
k2k
.
Baza interpolării este:
ncos,nsin,,cos,sin,2
1 .
Polinomul trigonometric de interpolare este:
ncosansinbcosasinb2
aP nn11
0
n
Condiţiile de interpolare conduc la un sistem cu matrice ortogonală cu soluţia:
n2
0k
01n2
k2f
1n2
2a
n:1j,1n2
jk2sin
1n2
k2f
1n2
2jsinf
1n2
2b
n2
1k
n2
1k
kkj
.1n2
jk2cos
1n2
k2f
1n2
2jcosf
1n2
2a
n2
0k
n2
0k
kkj
2.Funcţii spline de interpolare
Interpolarea polinomială, deşi este simplă, nu este stabilă numeric, limitând drastic gradul
polinomului de interpolare. Aceasta determină o aproximare locală, pe intervale, a funcţiei, prin
polinoame de grad scăzut (de obicei 3). In noduri se impun condiţii de interpolare şi de
continuitate (pentru funcţie şi derivatele ei).
Astfel pe suportul x0, x1,…,xn se consideră intervalele
[x0,x1), [x1,x2), …, [xn-1,xn)
Un spline cubic de interpolare este definit ca:
Si: [xi, xi+1) R, Si(x)=ai+bi(x-xi)+ci(x-xi)2+di(x-xi)
3
Page 5
Sistemul de coordonate local poate fi parametrizat prin schimbarea de variabilă:
i
i
i1i
i
h
xx
xx
xxt
Splineul în baza (1 t t2 t3) este:
si: [0,1) R, si(t)=i+it+it2+it
3
Utilizarea bazei Bernstein (1-t)3, 3t(1-t)2, 3t2(1-t), t3 simplifică în mod
considerabil calculul parametrilor splineurilor:
si: [0,1) R, si(t)=ai(1-t)3+3bit(1-t)
2+3cit2(1-t)+dit
3
Splineurile în clasă C1 utilizează baza Bernstein şi impun condiţii de interpolare de tip
Hermite:
Si(xi)=f(xi), i=0:n-1
Si’(xi)=f’(xi)
Sn-1(xn)=f(xn)
S’n-1(xn)=f’(xn)
şi condiţii de continuitate şi derivabilitate în nodurile interne:
Si(xi+1)=Si+1(xi+1), i=0:n-2
Si’(xi+1)=S’i+1(xi+1)
Parametrii funcţiei spline Si sunt:
1ii
1ii
1ii
ii
ii
ii
xfd
xf3
hxfc
,xf3
hxfb
,xfa
Splineurile în clasă C2 folosesc condiţii de interpolare de tip Lagrange:
Si(xi)=f(xi), i=0:n-1
Sn-1(xn)=f(xn)
şi condiţii de continuitate, derivabilitate şi curbură în nodurile interne:
Si(xi+1)=Si+1(xi+1), i=0:n-2
Si’(xi+1)=S’i+1(xi+1)
Si”(xi+1)=S”i+1(xi+1)
Pentru determinarea parametrilor mai sunt necesare două condiţii, care se pun la capete şi
sunt fie:
S0”(x0)=S”n-1(xn)=0 pentru splineuri naturale
S0’(x0)=f’(x0), S’n-1(xn)=f’(xn) pentru splineuri tensionate
Determinarea parametrilor funcţiilor spline impune rezolvarea unui sistem de ecuaţii cu
matrice tridiagonală.
Sistemul tridiagonal pentru funcţii spline naturale are forma:
Page 6
0
33
330
10000
200
002
0001
2
21
1
1
0
01
1
12
1
1
0
1122
1100
n
nn
n
nn
n
nnnnnh
aa
h
aa
h
aa
h
aa
c
c
c
c
hhhh
hhhh
ai=f(xi), i=0:n
1:1,23
11
nicc
h
h
aab ii
i
i
iii
1:0,3
1
nih
ccd
i
iii
Pentru funcţii spline tensionate avem sistemul:
1
1
0
21
1
1
0
01
1
12
0
0
01
1
1
0
11
1122
1100
00
33
33
33
33
2000
200
002
002
n
nn
n
nnnn
n
n
nn
nnnn
h
aaxf
h
aa
h
aa
h
aa
h
aa
xfh
aa
c
c
c
c
hh
hhhh
hhhh
hh
Probleme rezolvate 1. Se consideră funcția cunoscută prin tabelul:
x -1 -0.5 0.5 1
f(x) 1 0 2 1
a) Scrieţi expresia polinomul Lagrange de interpolare.
b) Scrieţi expresia polinomului de interpolare Newton şi expresia erorii.
c) Scrieţi o funcţie MATLAB, function c = coefNewton(x,y), care calculează coeficienţii
polinomului Newton de interpolare asociat funcţiei f cunoscută prin valorile y de coordonate x.
Soluție:
a)Expresia polinomului Lagrange de interpolare:
Page 7
5.1
)5.0)(5.0)(1(
75.0
)1)(5.0)(1(2
5.1
)1)(5.0)(5.0(
)5.01)(5.01)(11(
)5.0)(5.0)(1()1(
)15.0)(5.05.0)(15.0(
)1)(5.0)(1()5.0(
)15.0)(5.05.0)(15.0(
)1)(5.0)(1()5.0(
)11)(5.01)(5.01(
)1)(5.0)(5.0()1(
,0
,0
0
xxxxxxxxx
xxxf
xxxf
xxxf
xxxfn
kiiix
kx
n
kiiixx
n
k kxfxnP
b)Expresia polinomului de interpolare Newton
Pn(x)=f(x0)+(x-x0)F1[x0,x1]+(x-x0)(x-x1)F2[x0,x1,x2]+…
+(x-x0)…(x-xn-1)Fn[x0,x1,…,xn](n=3)
P3(x) = f(-1) + (x+1)F1[-1,-0.5]+ (x+1)(x+0.5)F2[-1,-0.5,0.5]
+(x+1)(x+0.5)(x-0.5)F3[-1,-0.5,0.5,1]
xi F0[xi] F1[xi,xi+1] F2[xi,xi+1,xi+2] F3[xi,xi+1,xi+2,xi+3]
-1 1 2
5.01
)5.0()1(
ff
3
8
5.01
22
3
8
11
3
8
3
8
-0.5 0
0.5 2 2
5.05.0
)5.0()5.0(
ff
1 1 2
15.0
)1()5.0(
ff
3
8
15.0
22
P3(x) = 1 + (x+1)(-2)+ (x+1)(x+0.5)3
8+(x+1)(x+0.5)(x-0.5) )
3
8( .
Eroarea interpolării are expresia:
!4
|)(||1)-0.5)(x-0.5)(x1)(x(x|]x,x,x,x[x,)Fx-(x)x-(x
)4(
32104n03
fxE ,
[-1,1].
c) Pn(x)=f(x0)+(x-x0)F1[x0,x1]+(x-x0)(x-x1)F2[x0,x1,x2]+…
+(x-x0)…(x-xn-1)Fn[x0,x1,…,xn] = c0+ c1 x1+ c2 x
2+...+ cn x
n.
function c = coefNewton(x,y)
z=difDiv(x,y);
n=length(x);
A=zeros(n,n);
A(1,n)=z(1);
c=zeros(n,1);
for i=2:n
A(i,n-i+1:n)=poly(x(1:i-1))*z(i);
Page 8
endfor
for i=1:n
c(i) = sum(A(:,i));
endfor
endfunction
function z=difDiv(x,y)
x=x(:);
y=y(:);
z=y;
n=length(y);
for i=1:n
z(i+1:n)=(z(i+1:n)-z(i))./(x(i+1:n)-x(i));
endfor
endfunction
2. Se consideră funcția cunoscută prin tabelul:
x -3 -1 0 4
f 2 0 3 6
a) Determinaţi expresia polinomului Lagrange de interpolare
b) Determinaţi expresia polinomului Newton de interpolare
c) Eroarea polinomului Newton de interpolare
a)
;))1(4)(04))(3(4(
))1()(0))(3((6
)40))(1(0))(3(0(
)4))(1())(3((3
)41)(01))(3(1(
)4)(0))(3((0
)43)(03))(1(3(
)4)(0))(1((2
xxx
xxxxxxxxxP
b)
;]],...,,[[],...,,[
],...,,[
;][][
],[
);(][
0
211101
10
10
1000101
000
p
pppp
ppxx
xxxFxxxFxxxF
xx
xFxFxxF
xfxF
X F0 F1 F2 F3
-3
-1
0
4
2
0
3
6
-1
3
¾
4/3
-9/20
-127/420
Page 9
F1: (2-0)/(-3+1)=-1; F2: (-1-3)/(-3-0)=4/3;
(0-3)/(-1-0)=3; (3-3/4)/(-1-4)=-9/20;
(3-6)/(0-4)=3/4 F3: (4/3+9/20)/(-3-4)=-107/420;
P=2+ (-1)(x-(-3))+4/3(x+3)(x+1) + (-127/420) (x+3)(x+1)x;
c)Eroarea
;4
)()4()1)(3()(
;)1(
)())...()(()(
)4(
)1(
10
fxxxxxE
n
fxxxxxxxE
n
nn
3. Pentru funcţia R]3,5[:f , 1x2)x(f , determinaţi polinomul generalizat de
interpolare Cebâşev de grad 2.
Soluţie. Facem schimbarea de variabilă
4
1xt1t4x
R)t(F]1,1[t:F
)t(Ta)t(Ta2
a)t(P 2211
0
2
Condiţiile de interpolare sunt:
)()()()(
)()()()(
)()()()(
2222211022
1122111012
0022011002
tFtTatTaatP
tFtTatTaatP
tFtTatTaatP
unde
2:0 ,6
12cos
i
iti .
2
3
6
5cos
02
cos
2
3
6cos
2
1
0
t
t
t
1323
5cos
6
5cos
1cos2
cos
1323
cos6
cos
210
210
210
aaa
aaa
aaa
a0=-1, a1=4 , a2=0 )(42
1)( 12 tTtP
3. Fie f o funcţie periodică impară, cu perioada 2, cunoscută prin valorile:
X /4 2/4 3/4
F(x) -1 0 1
a) Determinaţi polinomul de interpolare trigonometrică de grad 3: P3=b1sinx + b2sin2x + b3sin3x
Soluţie. a) Pentru a determina polinomul P3 este necesar să determinăm coeficienţii b1, b2, b3.
Page 10
Pentru aceasta avem:
14
3sin
4
2sin
4sin
43213
bbbP
04
6sin
4
4sin
4
2sin
4
23213
bbbP
14
9sin
4
6sin
4
3sin
4
33213
bbbP
De aici obţinem sistemul:
2222
0
2222
321
31
321
bbb
bb
bbb
care are soluţiile 0b1 , 1b2 şi 0b3 . Deci polinomul de interpolare trigonometric
este P3=-sin(2x).
4. Fie funcţia cunoscută prin tabelul:
X -1 0 1 2
f(x) -2 -1 2 4
a) Calculaţi funcţia spline cubică de clasă C1 care interpolează funcţia dată.
b) Calculaţi funcţia spline cubică de clasă C2 (spline natural) care interpolează funcţia dată.
Soluţie. a) spline în clasă C1:
Considerăm funcţii de interpolare liniară, locale pe subintervalele: [-1, 0], [0, 1],
[1, 2], de formaiii bxa)x(p
- condiţii de interpolare:
)()(
)()(
1 nnn
iii
xfxp
xfxp
4)2(
2)1(
1)0(
2)1(
2
2
1
0
p
p
p
p
- condiţii de racordare: )()( 111 iiii xpxp
2)1(
1)0(
1
0
p
p
2
1
42
2
1
2
11
0
22
22
1
00
ba
b
ba
ba
b
ba
0
2
1
3
1
1
2
2
1
1
0
0
b
a
b
a
b
a
Page 11
b) spline în clasă C2:
)(62
)(3)(2
)()()(
''
2'
32
iiii
iiiiii
iiiiiiii
xxdcs
xxdxxcbs
xxdxxcxxbas
3
2
2
2222
3
1
2
1111
3
0
2
0000
)1()1()1()(
)(
)1()1()1()(
xdxcxbaxs
xdxcxbaxs
xdxcxbaxs
- condiţii de interpolare:
4)2(
2)1(
1)0(
2)1(
2
2
1
0
s
s
s
s
- condiţii de racordare:
)1()1(
)0()0(
)1()1(
)0()0(
)1()1(
)0()0(
''
2
''
1
''
1
''
0
'
2
'
1
'
1
'
0
21
10
ss
ss
ss
ss
ss
ss
- condiţii pentru funcţii spline naturale:
0)2(
0)1(
''
2
''
0
s
s
02
02
262
262
32
32
4
2
1
2
2
0
211
100
2111
1000
21111
10000
2222
2
1
0
c
c
cdc
cdc
bdcb
bdcb
adcba
adcba
dcba
a
a
a
03
03
32
3
3
1
2
0
0
2
1
2
11
01
2111
100
111
00
22
2
0
2
1
0
dc
dc
bdcb
bdb
dcb
db
db
c
c
a
a
a
2
9
2
1
2
1
2
1
2
2
1
0
2
1
0
d
d
d
a
a
a
;
2
3
0
0
2
13
2
2
1
1
2
0
2
1
0
c
c
c
b
b
b
Page 12
Probleme Propuse 1. Calculaţi funcţiile spline în clasa C1 tensionate pentru funcţia cunoscută prin: x 1 2 3
F -1 1 18
f’
-3 9 27
2. Pentru nodurile A(-1,3), B(1,1),C(3,-3), D(5,0) calculaţi polinomul de interpolare sub formă
Lagrange şi Newton şi daţi expresia erorii interpolării.
3. Se consideră nodurile A(-1,1),B(2,3) şi C(4,5) prin care trec două splineuri naturale de
interpolare în clasă C2. Determinaţi-le.
4. Se consideră nodurile A(-1,1), B(2,3) şi C(4,5) prin care trec două splineuri de interpolare în
clasă C1 şi au pantele m(A)=1, m(B)=0, m(C)=-1. Determinaţi-le.
5. O funcţie este cunoscută prin tabloul:
x x(1) x(2) x(3)
y y(1) y(2) y(3) Scrieţi o funcţie Matlab pentru calculul funcţiilor spline cubice s1(x) şi s2(x)în clasă C2, naturale
cu semnătura function a=splC2n(x,y), în care a este o matrice cu 2 linii şi 4 coloane, conţinând
cei 4 parametri ai celor două splineuri.
6.Pentru funcţia )(xf cunoscută prin tabelul următor, calculaţi funcţiile spline )(0 xs şi
)(1 xs ,tensionate de clasă C2, dacă 2)1(''0 s şi .1)3(''1 s
x 1 2 3
f(x) 3 4 2
7. Calculaţi funcţiile spline )(0 xs şi )(1 xs în clasă C1 pentru funcţia )(xf cunoscută prin
tabelul următor:
x 1 2 4
f(x) 3 4 6
f’(x) 0 2 5
8. Calculați funcțiile spline C1 și spline C
2 naturale pentru funcția cunoscută prin :
x 1 2 3
f 5 -1 4
f’ 1 0 -1
