Metody Numeryczne Wykład 4 Wykład 5
Interpolacja wielomianowa Sformułowanie zadania interpolacji
Niech D ⊂ R i niech F b¸edzie pewnym zbiorem funkcji f : D → R. Niech x 0 , x 1 , . . . , x n
b¸edzie ustalonym zbiorem parami różnych punktów z D, zwanych w¸ ezłami.
Mówimy, że wielomian w interpoluje funkcj¸e f ∈ F w w¸ezłach x j , gdy w(x j ) = f (x j ), 0 ≤ j ≤ n.
Oznaczmy przez Π n przestrzeń liniow¸ a wielomianów stopnia co najwyżej n o współczyn- nikach rzeczywistych,
Π n = {w(x) = a n x n + a n−1 x n−1 + . . . + a 1 x + a 0 : a j ∈ R, 0 ≤ j ≤ n}.
Twierdzenie
Dla dowolnej funkcji f : D → R, istnieje dokładnie jeden wielomian interpolacyjny w ∈ Π n interpoluj¸ acy f w w¸ezłach x j , 0 ≤ j ≤ n.
Dowód
Wybierzmy w Π n dowoln¸ a baz¸e wielomianów φ j , 0 ≤ j ≤ n, Wtedy
Π n = span{φ 0 , φ 1 , . . . , φ n }.
gdzie span oznacza powłok¸e liniow¸ a, czyli zbiór wszystkich kombinacji liniowych funkcji φ j , 0 ≤ j ≤ n.
Każdy wielomian z Π n można wi¸ec jednoznacznie przedstawić w postaci rozwini¸ecia wzgl¸edem wybranej bazy. Warunkiem koniecznym i dostatecznym na to, aby wielomian w(·) = P n
j=0 c j φ j (·) interpolował f, jest spełnienie układu (n + 1) równań liniowych
n
X
j=0
c j φ j (x i ) = f (x i ), 0 ≤ i ≤ n.
z n + 1 niewiadomymi c j , który w postaci macierzowej ma postać
φ 0 (x 0 ) φ 1 (x 0 ) . . . φ n (x 0 ) φ 0 (x 1 ) φ 1 (x 1 ) . . . φ n (x 1 )
.. .
φ 0 (x n ) φ 1 (x n ) . . . φ n (x n )
c 0 c 1 .. . c n
=
f (x 0 ) f (x 1 )
.. . f (x n )
.
Aby wykazać, że ten układ ma jednoznaczne rozwi¸ azanie wystarczy, aby wektor zerowy był jedynym rozwi¸ azniem układu jednorodnego. Zauważmy, że układ jednorodny od- powiada interpolacji danych zerowych V i f (x i ) = 0. Istnienie niezerowego rozwi¸ azania byłoby wi¸ec równoważne istnieniu niezerowego wielomianu stopnia nie wi¸ekszego od n, który miałby n + 1 różnych zer x i , co jest niemożliwe, a wi¸ec wektor zerowy jest jedynym rozwi¸ azaniem układu jednorodnego, co mieliśmy wykazać.
W danym modelu obliczeniowym możemy wi¸ec wyznaczać jednoznacznie współczynniki c j wielomianu w.
Dokładniej Jeżeli (φ j ) n j=0 jest baz¸ a przestrzeni Π n wielomianów stopnia co najwyżej n, to zadanie iterpolacji wielomianowej polega na obliczeniu dla danej funkcji f współczynni- ków c j takich, że wielomian
w(·) =
n
X
j=0
c j φ j (·)
interpoluje f w punktach x j , 0 ≤ j ≤ n.
Wybór bazy wielomianowej
Można powiedzieć, że trudność zadania interpolacji jest zdeterminowana przez rozwi¸ azanie powyższego układu i zależy w istotny sposób od wybranej bazy (φ j ) n j=0 . W naturalny wi¸ec sposób powstaje problem wyboru ”wygodnej”bazy w Π n . Rozpatrzymy trzy bazy:
Lagrange, pot¸egow¸ a i Newtona.
BAZA LAGRANGE (KANONICZNA) Zdefiniujemy dla 0 ≤ j ≤ n wielomiany
l j (x) = (x − x 0 )(x − x 1 ) . . . (x − x j−1 )(x − x j+1 ) . . . (x − x n ) (x j − x 0 )(x j − x 1 ) . . . (x j − x j−1 )(x j − x j+1 ) . . . (x j − x n ) . Zauważmy że każdy wielomian l j jest stopnia dokładnie n oraz
l j (x i ) = 0 i 6= j, 1 i = j.
St¸ ad widać, że wielomiany te stanowi¸ a baz¸e w Π n , któr¸ a nazywamy baz¸ a Lagrange.
Wielomian interpolacyjny dla funkcji f zwany wielomianem interpolacyjnym Lagrange,
można w tym przypadku zapisać jako w(·) =
n
X
j=0
f (x j )l j (·).
Przykład
Skonstruujemy wielomian interpolacyjny Lagrange dla danych (1, 1), (2, 3), (4, 3).
l 0 (x) = a(x − 2)(x − 4), bo musi si¸e zerować w punktach x = 2, x = 4.
Ż¸ adaj¸ ac l 0 (1) = 1 - otrzymujemy równanie na nieznany współczynnik a, a(1 − 2)(1 − 4) = 1, a = 1 3 i l 0 (x) = 1 3 (x − 2)(x − 4). Podobnie określamy
l 1 (x) = − 1
2 (x − 1)(x − 4) ,
l 2 (x) = 1
6 (x − 1)(x − 2).
St¸ ad otrzymujemy nast¸epuj¸ ac¸ a postać wielomianu Lagrange:
w 2 (x) = f (x 0 )
3 (x − 2)(x − 4) − f (x 1 )
2 (x − 1)(x − 4) + f (x 2 )
6 (x − 1)(x − 2)
= 1
3 (x − 2)(x − 4) − 3
2 (x − 1)(x − 4) + 3
6 (x − 1)(x − 2)
Do obliczenia współczynników wielomianu Lagrange i jego wartości w zadanym punk- cie p, możemy napisać skrypt funkcyjny w OCTAVE na przykład o nazwie lagrange.m f unction lagrange(x, y, p)
n = length(x);
w = 0;
f or k = 1 : n b(k) = 1;
d(k) = 1;
f or j = 1 : n if j ∼ k
b(k) = b(k) ∗ (x(k) − x(j));
d(k) = d(k) ∗ (p − x(j));
end end
c(k) = y(k)/b(k);
w = w + c(k) ∗ d(k);
end
c;w
BAZA POTE ¸ GOWA (NATURALNA)
Znacznie prościej można obliczyć wartość wielomianu interpolacyjnego (a także jego po- chodnych), gdy jest on dany w najcz¸eściej używanej bazie pot¸egowej, V
φ j (x) = x j . Jeśli
w(x) = a 0 + a 1 x + . . . + a n x n ,
to również
w(x) = (. . . (a n x + a n−1 )x + a n−2 )x + . . . + a 1 )x + a 0 ,
co sugeruje zastosowanie nast¸epuj¸ acego schematu Hornera do obliczenia wartości wielo- mianu w(x).:
v n = a n ;
f or j = (n − 1) : 0 v j = v j+1 ∗ x + a j , end
Po wykonaniu tego algorytmu w(x) = v 0 .
Schemat Hornera wymaga wykonania tylko n mnożeń i n dodawań.
Zauważmy jednak ˙że w przypadku bazy pot¸egowej macierz (x j i ) n i,j=0 układu jest pełna.
Jest to tzw. macierz Vandermonde. Obliczenie współczynników wielomianu interpolacyj- nego w bazie pot¸egowej bezpośrednio z tego układu jest kosztowne bo wymaga n 3 operacji.
Przykład
Znajdziemy wielomian interpolacyjny w bazie pot¸egowej w 2 (x) = c 0 + c 1 x + c 2 x 2 ”prze- chodz¸ acy”przez punkty (1, 1), (2, 3) i (4, 3).
Warunki interpolacji
w 2 (x 0 ) = c 0 + 1c 1 + 1c 2 = 1, w 2 (x 1 ) = c 0 + 2c 1 + 4c 2 = 3, w 2 (x 2 ) = c 0 + 4c 1 + 16c 2 = 3.
Postać macierzowa tego układu równań liniowych
1 1 1 1 2 4 1 4 16
c 0 c 1 c 2
=
1 3 3
Rozwi¸ azuj¸ ac ten układ za pomoc¸ a OCTAVE
A = [1 1 1; 1 2 4; 1 4 16];
b = [1; 3; 3];
c = A\b
Uwzgl¸edniaj¸ ac bł¸edy zaokr¸ agleń otrzymujemy c 0 = − 7 3 , c 1 = 4, c 2 = − 2 3 . Wielomian interpolacyjny ma wi¸ec postać
w 2 (x) = (−2x 2 + 12x − 7)/3
BAZA NEWTONA
Rozwi¸ azaniem pośrednim, które ł¸ aczy prostot¸e obliczenia współczynników z prostot¸ a ob- liczenia wielomianu w(x) jest wybór bazy Newtona.
p 0 (x) = 1,
p j (x) = (x − x 0 )(x − x 1 ) . . . (x − x j−1 ), 1 ≤ j ≤ n.
W tym przypadku współczynniki rozwini¸ecia wielomianu interpolacyjnego b¸edziemy ozna- czać przez b j , zaś wielomian Newtona ma postać
w =
n
X
j=0
b j p j
Zwróćmy uwag¸e na ważn¸ a własność wielomianu Newtona, wynikaj¸ a z konstrukcji jego bazy.
Jeśli w j ∈ Π j jest wielomianem interpolacyjnym dla funkcji f opartym na w¸ezłach x 0 , x 1 , . . . , x j , 0 ≤ j ≤ n,to
w 0 = b 0
w j = w j−1 + b j p j , 1 ≤ j ≤ n
Wartość w(x) można obliczać, stosuj¸ ac prost¸ a modyfikacj¸e algorytmu Hornera:
v n = b n ;
f orj = n − 1 : 0;
v j = v j−1 ∗ (x − x j ) + b j . end
Ponadto układ równań jest trójk¸ atny dolny ( tzn.p j (x i ) = 0 dla i < j),co pozwala na obliczenie algorytmem, który Państwo poznaliście na wykładzie 3 o koszcie n 2 ,
b 0 = f (x 0 );
f or j = 1 : n
b j =
f (x j ) − P j−1
s=0 b s p s (x j )
/p j (x j ).
Zwykle jednak do obliczania współczynników b j stosuje inny algorytm, który teraz przed- stawimy.
Algorytm różnic dzielonych
Różnic¸e dzielon¸ a funkcji f opart¸ a na różnych w¸ azłach t 0 , t 1 , . . . , t s , gdzie s ≥ 1, definiuje si¸e indukcyjnie jako
f [t 0 , t 1 , . . . , t s ] = f [t 1 , t 2 , . . . , t s ] − f [t 0 , t 1 , . . . , t s−1 ] t s − t 0
Zachodzi nast¸epuj¸ ace ważne twierdzenie o współczynnikach wielomianu interpolacyjnego Newtona.
Twierdzenie
Współczynniki wielomianu interpolacyjnego Newtona dla danej funkcji f dane s¸ a przez różnice dzielone f w w¸ezłach x 0 , x 1 , . . . , x j , tzn.
b j = f [x 0 , x 1 , . . . , x j ], 0 ≤ j ≤ n Dowód
Dla 0 ≤ i ≤ j ≤ n, oznaczmy przez w i,j wielomian z Π j−i interpoluj¸ acy f w w¸ezłach x i , x i+1 , . . . , x j . Wtedy ma miejsce nast¸epuj¸ aca równość dla i < j:
^ x w i,j (x) = (x − x i )w i+1,j − (x − x j )w i,j−1 (x) (x j − x i ) ,
Aby wykazać t¸e równość wystarczy pokazać, że prawa jej strona, któr¸ a oznaczymy przez v(x), przyjmuje wartości f (x s ) dla x = x s , i ≤ s ≤ j.
Jeśli i + 1 ≤ s ≤ j − 1 to
v(x s ) = (x s − x i )f (x s ) − (x s − x j )f (x s )
x j − x i = f (x s ).
Ponadto
v(x i ) = −(x i − x j )
x j − x i f (x i ) = f (x i ) oraz
v(x j ) = (x j − x i ) x j − x i
f (x j ) = f (x j )
Wykazaliśmy, że v jest wielomianem z Π j−i interpoluj¸ acym f w w¸ezłach x s , i ≤ s ≤ j, czyli w i,j = v.
Dalej post¸epujemy indukcyjnie ze wzgl¸edu na stopień wielomianu interpolacyjnego.
Dla n = 0 mamy b 0 = f (x 0 . Niech n ≥ 1. Jak wiemy z konstrukcji wielomianu Newtona w 0,n = w 0,n−1 (x) + b n p n (x)
i z założenia indukcyjnego b j = f [x 0 , x 1 , . . . , x j ], dla 0 ≤ j ≤ n − 1
Aby pokazać podobn¸ a równość dla b n wykorzysamy równość na w i,j z i = 0, j = n, która ma wtedy postać
w 0,n (x) = (x − x 0 )w 1,n (x) − (x − x n )w 0,n−1 (x) x n − x 0
Zauważmy, że b n jest współczynnikiem przy x n w wielomianie w 0,n . Z zal ożenia induk- cyjnego wynika, że współczynniki przy x n−1 w wielomianach w 1,n , w 0,n−1 s¸ a ilorazami różnicowymi opartymi odpowiednio na w¸ezłach x 1 , x 2 , . . . , x n ,i x 0 , x 1 , . . . x n−1 .
St¸ ad
b n = f [x 1 , x 2 , . . . , x n ] − f [x 0 , x 1 , . . . , x n−1 ] x n − x 0
= f [x 0 , x 1 , . . . , x n ] co mieliśmy udowodnić.
Różnic¸e dzielon¸ a można obliczć na podstawie wartości f (x j ), 0 ≤ j ≤ n, buduj¸ ac tabelk¸e x 0 f [x 0 ]
x 1 f [x 1 ] f [x 0 , x 1 ] x 2 f [x 2 ] f [x 1 , x 2 ] .. .
x n f [x n ] f [x n−1 , x n ] . . . f [x 0 , x 1 , . . . , x n ].
Zauważmy, że po drodze obliczamy f [x i , x i+1 , . . . , x j ] dla wszystkich 0 ≤ i < j ≤ n, a wi¸ec w szczególności interesuj¸ ace nas różnice dzielone f [x 0 , x 1 , . . . , x n ].
Przykład
Znajdziemy wielomian Newtona interpoluj¸ acy dane: (1, 1), (2, 3), (4, 3).
Obliczamy różnice dzielone f [x 0 , x 1 ] = f [x x1]−f [x
0]
1
−x
0= 3−1 2−1 = 2 f [x 1 , x 2 ] = f [x x2]−f [x
1]
2
−x
1= 3−3 4−2 = 0
f [x 0 , x 1 , x 2 ] = f [x1,x x
2]−f [x
1,x
0]
2