MN 09 Interpol i Aproks, metody numeryczne


MN09

Interpolacja wielomianowa

Zadanie interpolacji, czyli poprowadzenia krzywej zadanego rodzaju przez zestaw danych punktów, jest jednym z podstawowych zadań obliczeniowych. Stosuje się je nagminnie w najróżniejszych dziedzinach życia, np. wtedy, gdy trzeba

Interpolację stosuje się szczególnie chętnie w samej numeryce. Na przykład idea metody siecznych polega na tym, by funkcję, której miejsca zerowego szukamy, przybliżyć prostą interpolującą tę funkcję w dwóch punktach. Metody numerycznego całkowania oraz rozwiązywania równań różniczkowych także korzystają z interpolacji.

Wielomian 0x01 graphic
(czerwony) stopnia 6, interpolujący 7 zadanych wartości (zaznaczone na zielono) danej funkcji 0x01 graphic
.

Niech 0x01 graphic
i niech 0x01 graphic
będzie pewnym zbiorem funkcji 0x01 graphic
. Niech 0x01 graphic
będzie ustalonym zbiorem parami różnych punktów z 0x01 graphic
, zwanych później węzłami. Powiemy, że wielomian 0x01 graphic
interpoluje funkcję 0x01 graphic
w węzłach 0x01 graphic
, gdy

0x01 graphic

Oznaczmy przez 0x01 graphic
przestrzeń liniową wielomianów stopnia co najwyżej 0x01 graphic
o współczynnikach rzeczywistych,

0x01 graphic

Zadanie znalezienia wielomianu interpolującego zadane wartości nazywamy zadaniem interpolacji Lagrange'a.

Twierdzenie o istnieniu i jednoznaczności zadania interpolacji Lagrange'a. Dla dowolnej funkcji 0x01 graphic
istnieje dokładnie jeden wielomian 0x01 graphic
interpolujący 0x01 graphic
w węzłach 0x01 graphic
, 0x01 graphic

Dowód. Wybierzmy w 0x01 graphic
dowolną bazę wielomianów 0x01 graphic
, 0x01 graphic
, tak, że 0x01 graphic
Wtedy każdy wielomian z 0x01 graphic
można jednoznacznie przedstawić w postaci rozwinięcia względem wybranej bazy. Warunkiem koniecznym i dostatecznym na to, aby wielomian

0x01 graphic

interpolował 0x01 graphic
jest spełnienie układu 0x01 graphic
równań liniowych

0x01 graphic

z 0x01 graphic
niewiadomymi 0x01 graphic
, który w postaci macierzowej wygląda następująco:

0x01 graphic

Aby wykazać, że układ ten ma jednoznaczne rozwiązanie wystarczy, aby wektor zerowy był jedynym rozwiązaniem układu jednorodnego. Rzeczywiście, układ jednorodny odpowiada interpolacji danych zerowych, 0x01 graphic
, 0x01 graphic
. Istnienie niezerowego rozwiązania byłoby więc równoważne istnieniu niezerowego wielomianu stopnia nie większego od 0x01 graphic
, który miałby 0x01 graphic
różnych zer 0x01 graphic
, co jest niemożliwe.

Zadanie znalezienia dla danej funkcji 0x01 graphic
jej wielomianu interpolacyjnego stopnia co najwyżej 0x01 graphic
jest więc dobrze zdefiniowane, tzn. rozwiązanie istnieje i jest wyznaczone jednoznacznie. Zauważmy, że wielomian interpolacyjny 0x01 graphic
jako taki nie może być wynikiem obliczeń w naszym modelu obliczeniowym. Możemy natomiast wyznaczyć jego współczynniki 0x01 graphic
w wybranej bazie.

Definicja. Niech 0x01 graphic
będzie bazą w przestrzeni 0x01 graphic
wielomianów stopnia co najwyżej 0x01 graphic
. Zadanie interpolacji wielomianowej polega na obliczeniu dla danej funkcji 0x01 graphic
współczynników 0x01 graphic
takich, że wielomian

0x01 graphic

interpoluje 0x01 graphic
w punktach 0x01 graphic
, 0x01 graphic
.

Wybór bazy wielomianowej

Jak już wiemy, zadanie interpolacji Lagrange'a sprowadza się do rozwiązania układu równań liniowych. Okazuje się, że w zależności od wyboru sposobu reprezentacji naszego wielomianu (czyli od wyboru bazy wielomianowej 0x01 graphic
), układ ten może być albo bardzo łatwy do rozwiązania, albo bardzo trudny. Co więcej, jego rozwiązanie w arytmetyce 0x01 graphic
może napotykać na większe bądź mniejsze trudności (w zależności np. od uwarunkowania macierzy układu, który musimy rozwiązać).

W naturalny sposób powstaje więc problem wyboru "wygodnej" bazy w 0x01 graphic
. Rozpatrzymy trzy bazy: Lagrange'a, potęgową i Newtona.

Baza Lagrange'a (kanoniczna)

Zdefiniujmy dla 0x01 graphic
wielomiany

0x01 graphic

Zauważmy, że każdy z 0x01 graphic
jest stopnia dokładnie 0x01 graphic
oraz

0x01 graphic

Teraz widać, że wielomiany te stanowią bazę w 0x01 graphic
, którą nazywamy bazą Lagrange'a. Macierz układu zadania interpolacji jest w takim wypadku identycznością i w konsekwencji 0x01 graphic
, 0x01 graphic
. Wielomian interpolacyjny dla funkcji 0x01 graphic
można więc zapisać jako

0x01 graphic

Koszt kombinatoryczny rozwiązania zadania interpolacji jest przy tym zerowy.

Wzory barycentryczne

Przypuśćmy, że chcielibyśmy obliczyć wartość wielomianu interpolacyjnego 0x01 graphic
w punkcie 0x01 graphic
różnym od 0x01 graphic
, 0x01 graphic
. Podstawiając

0x01 graphic

oraz 0x01 graphic
mamy pierwszy wzór barycentryczny

0x01 graphic

i ostatecznie dostajemy tzw. drugi wzór barycentryczny na wielomian interpolacyjny,

0x01 graphic

gdzie 0x01 graphic
. W ostatniej równości wykorzystaliśmy fakt, że 0x01 graphic
, co łatwo widzieć, rozpatrując zadanie interpolacji funkcji 0x01 graphic
. Drugi wzór barycentryczny jest korzystniejszy w implementacji.

Dla wielu układów węzłów, wagi 0x01 graphic
są zadane jawnymi wzorami, np. dla węzłów równoodległych (niezależnie od tego, na jakim odcinku) wagi w drugim wzorze barycentrycznym wynoszą po prostu

0x01 graphic

Również dla węzłów Czebyszewa istnieją eleganckie wzory na takie współczynnki.

Można pokazać, że wartość 0x01 graphic
wielomianu interpolacyjnego obliczona w arytmetyce 0x01 graphic
według pierwszego wzoru barycentrycznego spełnia

0x01 graphic

gdzie 0x01 graphic
, a więc jest to algorytm numerycznie poprawny.

Zachowanie drugiej postaci wzoru barycentrycznego w arytmetyce 0x01 graphic
jest nieco bardziej skomplikowane.

Baza potęgowa (naturalna)

Znacznie prościej można obliczyć wartość wielomianu interpolacyjnego, (a także jego pochodnych), gdy jest on dany w najczęściej używanej bazie potęgowej, 0x01 graphic
, 0x01 graphic
. Jeśli bowiem

0x01 graphic

to również

0x01 graphic

co sugeruje zastosowanie następującego schematu Hornera do obliczenia 0x01 graphic
:

Algorytm Algorytm Hornera

0x01 graphic

0x01 graphic

for (j=n-1; j >= 0 ; j--)

0x01 graphic
;

Po wykonaniu tego algorytmu 0x01 graphic
. Schemat Hornera wymaga wykonania tylko 0x01 graphic
mnożeń i 0x01 graphic
dodawań. Ma on również głębszy sens, bo jego produktem ubocznym mogą być także wartości pochodnych naszego wielomianu w 0x01 graphic
. Algorytm Hornera okazuje się optymalny. Każdy inny algorytm obliczający dokładną wartość wielomianu, gdy danymi są współczynniki wielomianu, wymaga wykonania co najmniej 0x01 graphic
mnożeń i 0x01 graphic
dodawań. Algorytm Hornera jest też numerycznie poprawny.

Zauważmy jednak, że w przypadku bazy potęgowej macierz 0x01 graphic
układu zadania interpolacji jest pełna. Jest to tzw. macierz Vandermonde'a. Obliczenie współczynników wielomianu interpolacyjnego w bazie potęgowej bezpośrednio z tego układu, stosując jedną ze znanych nam już metod, kosztowałoby rzędu 0x01 graphic
operacji arytmetycznych. Co gorsza, w często spotykanym przypadku, gdy węzły interpolacji są równoodległe, ta macierz jest bardzo źle uwarunkowana.

Baza Newtona

Rozwiązaniem pośrednim, które łączy prostotę obliczenia współczynników z prostotą obliczenia wartości 0x01 graphic
i ewentualnie jego pochodnych, jest wybór bazy Newtona,

0x01 graphic

W tym przypadku współczynniki rozwinięcia wielomianu interpolacyjnego będziemy oznaczać przez 0x01 graphic
,

0x01 graphic

Zwróćmy od razu uwagę na ważną własność bazy Newtona. Jeśli 0x01 graphic
jest wielomianem interpolacyjnym dla funkcji 0x01 graphic
opartym na węzłach 0x01 graphic
, 0x01 graphic
, to 0x01 graphic
oraz

0x01 graphic

Wartość 0x01 graphic
możemy obliczyć, stosując prostą modyfikację algorytmu Hornera:

Algorytm Algorytm Hornera dla bazy Newtona

0x01 graphic

0x01 graphic

for (j=n-1; j >= 0 ; j--)

0x01 graphic
;

Ponadto układ równań zadania interpolacji jest trójkątny dolny, o specyficznej strukturze, dzięki czemu można stworzyć elegancki algorytm, który teraz przedstawimy.

Algorytm różnic dzielonych

Różnicę dzieloną funkcji 0x01 graphic
opartą na różnych węzłach 0x01 graphic
, gdzie 0x01 graphic
, definiuje się indukcyjnie jako

0x01 graphic

Zachodzi następujące ważne twierdzenie.

Twierdzenie O różnicach dzielonych. Współczynniki 0x01 graphic
wielomianu interpolacyjnego Newtona dla danej funkcji 0x01 graphic
dane są przez różnice dzielone 0x01 graphic
w węzłach 0x01 graphic
, tzn.

0x01 graphic

Dowód. Dla 0x01 graphic
, oznaczmy przez 0x01 graphic
wielomian z 0x01 graphic
interpolujący 0x01 graphic
w węzłach 0x01 graphic
. Wtedy ma miejsce następująca równość (0x01 graphic
):

0x01 graphic

Aby ją pokazać wystarczy, że prawa strona tej równości, którą oznaczymy przez 0x01 graphic
, przyjmuje wartości 0x01 graphic
dla 0x01 graphic
, 0x01 graphic
. Rzeczywiście, jeśli 0x01 graphic
to

0x01 graphic

Ponadto

0x01 graphic

oraz podobnie 0x01 graphic
. Stąd 0x01 graphic
jest wielomianem z 0x01 graphic
interpolującym 0x01 graphic
w węzłach 0x01 graphic
, 0x01 graphic
, czyli 0x01 graphic
.

Dalej postępujemy indukcyjnie ze względu na stopień 0x01 graphic
wielomianu interpolacyjnego. Dla 0x01 graphic
mamy oczywiście 0x01 graphic
. Niech 0x01 graphic
. Ponieważ, jak łatwo zauważyć,

0x01 graphic

z założenia indukcyjnego mamy 0x01 graphic
dla 0x01 graphic
. Aby pokazać podobną równość dla 0x01 graphic
, zauważmy, że

0x01 graphic

Zauważmy teraz, że 0x01 graphic
jest współczynnikiem przy 0x01 graphic
w wielomianie 0x01 graphic
. Z założenia indukcyjnego wynika, że współczynniki przy 0x01 graphic
w wielomianach 0x01 graphic
i 0x01 graphic
są ilorazami różnicowymi opartymi odpowiednio na węzłach 0x01 graphic
i 0x01 graphic
. Stąd

0x01 graphic

co kończy dowód.

Różnicę dzieloną 0x01 graphic
można łatwo obliczyć na podstawie wartości 0x01 graphic
, 0x01 graphic
, budując następującą tabelkę:

0x01 graphic

Wyznaczenie wielomianu 0x01 graphic
interpolującego zestaw punktów 0x01 graphic
algorytmem różnic dzielonych.

Zauważmy przy tym, że "po drodze" obliczamy 0x01 graphic
dla wszystkich 0x01 graphic
, a więc w szczególności również interesujące nas różnice dzielone 0x01 graphic
. Stąd i z twierdzenia o różnicach dzielonych wynika algorytm obliczania współczynników 0x01 graphic
wielomianu interpolacyjnego w bazie Newtona. Po wykonaniu następującego algorytmu,

Algorytm Metoda różnic dzielonych

0x01 graphic

for (j = 0; j <= n; j++)

0x01 graphic
= 0x01 graphic
;

for (j = 0; j <= n; j++)

for (k = n; k >= j; k--)

0x01 graphic
= 0x01 graphic
;

współczynniki 0x01 graphic
na końcu algorytmu zawierają współczynniki wielomianu interpolacyjnego w bazie Newtona.

Wyznaczenie tego samego wielomianu 0x01 graphic
, interpolującego zestaw punktów 0x01 graphic
algorytmem różnic dzielonych --- wykonanym tym razem in situ. Okazuje się, że przy realizacji w 0x01 graphic
algorytmu różnic dzielonych istotną rolę odgrywa porządek węzłów. Można pokazać, że --- o ile węzły są uporządkowane nierosnąco lub niemalejąco --- algorytm liczenia 0x01 graphic
jest numerycznie poprawny ze względu na dane interpolacyjne 0x01 graphic
, a cały algorytm różnic dzielonych daje w arytmetyce 0x01 graphic
współczynniki wielomianu interpolacyjnego, będące niewielkim zaburzeniem wartości dokładnych.

Uwarunkowanie

Danymi w zadaniu interpolacji są zarówno wartości interpolowanej funkcji, jak i węzły interpolacji. Traktując węzły jako sztywno zadane parametry zadania i dopuszczając jedynie zaburzenia wartości funkcji, można pokazać, że jeśli zamiast 0x01 graphic
rozpatrzyć jej zaburzenie 0x01 graphic
, gdzie 0x01 graphic
, to

0x01 graphic

gdzie

0x01 graphic

Znacznie rzadziej rozważa się uwarunkowanie zadania interpolacji ze względu na zaburzenie węzłów. Warto zaznaczyć, że zaburzenie danych interpolacji tylko w jednym punkcie może mieć wpływ na przebieg całego wielomianu interpolacyjnego, co ukazuje poniższy przykład:

Przykład. Pokażemy zmianę kilku bazowych wielomianów Lagrange'a stopnia 10 (dla węzłów równoodległych w 0x01 graphic
) w sytuacji, gdy trzeci węzeł interpolacji zostanie zaburzony o 0.01.

Wybrane wielomiany bazowe Lagrange'a oparte na węzłach równoodległych (zielone) kontra te same wielomiany, oparte na tych samych węzłach, z jednym wyjątkiem: węzeł 0x01 graphic
został zmieniony na 0x01 graphic
(czerwone).

Jak widać, to lokalne zaburzenie danych może powodować wyraźne globalne zaburzenie całego wielomianu interpolacyjnego.

MATLAB i Octave mają wbudowaną funkcję wyznaczającą wielomian, interpolujący zadane wartości: jeśli x jest wektorem zawierającym 0x01 graphic
węzłów, a y --- wektorem zawierającym wartości w węzłach, to

c = polyfit(x,y,N-1);

daje współczynniki wielomianu interpolacyjnego (Ostatni argument jest równy 0x01 graphic
, bo taki powinien być stopień wielomianu interpolacyjnego Lagrange'a).

Co ciekawe (i budzące trochę zgrozy) --- wielomian (zarówno w MATLABie, jak w Octave) jest wyznaczany w bazie naturalnej, przez rozwiązanie układu równań z macierzą Vandermonde'a, a więc w sposób najgorszy z możliwych. Napisz odpowiedni kod i wyślij do Octave-forge.

Aby teraz wyznaczyć wartości takiego wielomianu w zadanych punktach 0x01 graphic
, także musimy użyć specjalnej funkcji,

Y = polyval(c,X);

Implementuje ona algorytm Hornera.

Przykład. Interpolujemy tabelkę

0x01 graphic

2

1

0

0x01 graphic

5

2

1

wielomianem stopnia co najwyżej 2.

octave:1> x = [2, 1, 0]

x =

2 1 0

octave:2> y = [5, 2, 1]

y =

5 2 1

octave:3> c = polyfit(x,y,2)

c =

1 0 1

octave:4> polyval(c,3)

ans = 10

Zgodnie z przewidywaniami, otrzymaliśmy wielomian 0x01 graphic
. Wartość tego wielomianu dla 0x01 graphic
rzeczywiście jest równa 10.

A co się stanie, gdy będziemy szukać wielomianu stopnia niższego?

octave:6> c1 = polyfit(x,y,1)

c1 =

2.00000 0.66667

Też "coś" zostało obliczone --- wielomian 0x01 graphic
. Okazuje się, że to wielomian najlepiej pasujący do danych w sensie aproksymacji średniokwadratowej.

Warto jeszcze może wiedzieć, że polyfit można także wywołać dla jeszcze wyższego stopnia wielomianu, jednak, co niespodziewane, wynikiem nie będzie wielomian stopnia 2, uzyskany poprzednio:

octave:7> c3 = polyfit(x,y,3)

c3 =

0.21429 0.35714 0.42857 1.00000

Wynika to stąd, że gdy dopuszczalny stopień wielomianu jest wyższy niż wymagany w zadaniu interpolacji Lagrange'a, zadanie interpolacji ma nieskończenie wiele rozwiązań. Funkcja polyfit wybiera z nich to, które spełnia warunek, że norma euklidesowa wektora współczynników wielomianu jest najmniejsza z możliwych.

Przypadek węzłów wielokrotnych

Uogólnieniem rozpatrzonego zadania interpolacji jest zadanie interpolacji Hermite'a. Zakładamy, że oprócz (różnych) węzłów 0x01 graphic
dane są również ich krotności 0x01 graphic
, 0x01 graphic
, przy czym 0x01 graphic
. Należy skonstruować wielomian 0x01 graphic
taki, że

0x01 graphic

Oczywiście zakładamy przy tym, że odpowiednie pochodne funkcji 0x01 graphic
istnieją.

Lemat. Zadanie interpolacji Hermite'a ma jednoznaczne rozwiązanie.

Dowód. Istnienie i jednoznaczność rozwiązania można uzasadnić tak samo jak w przypadku węzłów jednokrotnych. Przedstawiając wielomian w dowolnej bazie otrzymujemy układ 0x01 graphic
równań z 0x01 graphic
niewiadomymi, który dla zerowej prawej strony ma jedynie rozwiązanie zerowe. Inaczej bowiem istniałby wielomian niezerowy stopnia nie większego niż 0x01 graphic
, który miałby zera o łącznej krotności większej niż 0x01 graphic
.

Nas oczywiście interesuje konstrukcja wielomianu 0x01 graphic
. W tym celu ustawimy węzły 0x01 graphic
w ciąg

0x01 graphic

i zdefiniujemy uogólnioną bazę Newtona w 0x01 graphic
jako

0x01 graphic

Uogólnimy również pojęcie różnicy dzielonej na węzły powtarzające się, kładąc

0x01 graphic

dla 0x01 graphic
, oraz

0x01 graphic

dla 0x01 graphic
. Zauważmy, że przy tej definicji różnice 0x01 graphic
możemy łatwo obliczyć stosując schemat podobny do tego z przypadku węzłów jednokrotnych.

Twierdzenie O różnicach dzielonych dla interpolacji Hermite'a

Współczynniki 0x01 graphic
wielomianu interpolacyjnego Hermite'a w bazie Newtona,

0x01 graphic

dane są przez odpowiednie różnice dzielone, tzn.

0x01 graphic

Dowód. Dowód przeprowadzimy podobnie jak dla węzłów jednokrotnych. Niech 0x01 graphic
oznacza wielomian interpolacyjny Hermite'a oparty na (być może powtarzających się) węzłach 0x01 graphic
. To znaczy, 0x01 graphic
interpoluje 0x01 graphic
w węzłach 0x01 graphic
takich, że 0x01 graphic
występuje w ciągu 0x01 graphic
, a jego krotność jest liczbą powtórzeń 0x01 graphic
w tym ciągu.

Zauważmy najpierw, że dla 0x01 graphic
zachodzi znany nam już wzór,

0x01 graphic

Rzeczywiście, oznaczmy przez 0x01 graphic
prawą stronę powyższej równości. Dla 0x01 graphic
mniejszego od krotności danego węzła 0x01 graphic
w ciągu 0x01 graphic
, mamy 0x01 graphic
, a ponieważ

0x01 graphic

to

0x01 graphic

Korzystając z tego wzoru sprawdzamy, że 0x01 graphic
spełnia odpowiednie warunki interpolacyjne, a stąd 0x01 graphic
.

Dalej postępujemy indukcyjnie ze względu na 0x01 graphic
. Dla 0x01 graphic
mamy 0x01 graphic
. Dla 0x01 graphic
wystarczy pokazać, że 0x01 graphic
. W tym celu rozpatrzymy dwa przypadki.

Jeśli 0x01 graphic
, to mamy jeden węzeł 0x01 graphic
o krotności 0x01 graphic
. Wielomian interpolacyjny jest wtedy postaci

0x01 graphic

a stąd 0x01 graphic
. Jeśli zaś 0x01 graphic
, to równość 0x01 graphic
wynika z wcześniej wyprowadzonych wzorów oraz z założenia indukcyjnego.

Zauważmy, ze pojęcie różnicy dzielonej formalnie zdefiniowaliśmy jedynie dla ciągu węzłów postaci 0x01 graphic
, gdzie 0x01 graphic
są parami różne. Tą definicję można rozszerzyć do dowolnego ciągu węzłów. Można bowiem powiedzieć, że 0x01 graphic
jest współczynnikiem przy 0x01 graphic
wielomianu 0x01 graphic
interpolującego 0x01 graphic
w węzłach 0x01 graphic
(uwzględniając krotności). Równoważnie,

0x01 graphic

Błąd interpolacji

Gdy mamy do czynienia z funkcją, która jest "skomplikowana", często dobrze jest zastąpić ją funkcją "prostszą". Mówimy wtedy o aproksymacji funkcji. Funkcję musimy również aproksymować wtedy, gdy nie jesteśmy w stanie uzyskać pełnej o niej informacji. Na przykład, gdy funkcja reprezentuje pewien proces fizyczny, często zdarza się, że dysponujemy jedynie ciągiem próbek, czyli wartościami tej funkcji w pewnych punktach. Jasne jest, że chcielibyśmy przy tym, aby błąd aproksymacji był możliwie mały.

Podobnie ma się sprawa w przypadku implementacji funkcji elementarnych (0x01 graphic
) w bibliotece funkcji matematycznych, czy wręcz w procesorze. Tam również najchętniej poszukiwalibyśmy sposobu taniego przybliżenia wartości dokładnej funkcji. I rzeczywiście, często w tym celu stosuje się m.in. specjalnie konstruowaną aproksymację wielomianową.

Z tego punktu widzenia, intepolacja wielomianowa może być traktowana jako jeden ze sposobów aproksymacji funkcji, opartym na próbkowaniu. Naturalnym staje się więc pytanie o błąd takiej aproksymacji.

Niech 0x01 graphic
będą (niekoniecznie różnymi) węzłami należącymi do pewnego (być może nieskończonego) przedziału 0x01 graphic
. Dla danej funkcji 0x01 graphic
, przez 0x01 graphic
rozważamy, tak jak w całym wykładzie, wielomian interpolacyjny stopnia co najwyżej 0x01 graphic
interpolujący 0x01 graphic
w zadanych węzłach. W przypadku węzłów wielokrotnych jest to oczywiście wielomian interpolacyjny Hermite'a; gdy węzły są jednokrotne, mamy do czynienia z interpolacją Lagrange'a.

Lemat Postać błędu interpolacji. Dla dowolnego punktu 0x01 graphic
błąd interpolacji w 0x01 graphic
wyraża się wzorem

0x01 graphic

Jeśli ponadto 0x01 graphic
, czyli pochodna 0x01 graphic
w 0x01 graphic
istnieje i jest ciągła, to

0x01 graphic

gdzie 0x01 graphic
jest pewnym punktem należącym do najmniejszego przedziału zawierającego punkty 0x01 graphic
.

Dowód. Możemy założyć, że 0x01 graphic
nie jest żadnym z węzłów 0x01 graphic
, 0x01 graphic
. Niech 0x01 graphic
będzie wielomianem interpolacyjnym funkcji 0x01 graphic
opartym na węzłach 0x01 graphic
i dodatkowo na węźle 0x01 graphic
. Mamy wtedy

0x01 graphic

a ponieważ z warunku interpolacyjnego 0x01 graphic
, to mamy też pierwszą równość w lemacie.

Aby pokazać drugą część lematu, rozpatrzmy funkcję 0x01 graphic
,

0x01 graphic

Z warunków interpolacyjnych na 0x01 graphic
wynika, że funkcja 0x01 graphic
ma punkty zerowe o łącznej krotności co najmniej 0x01 graphic
. Wykorzystując twierdzenie Rolle'a wnioskujemy stąd, że 0x01 graphic
ma zera o łącznej krotności co najmniej 0x01 graphic
, 0x01 graphic
ma zera o łącznej krotności co najmniej 0x01 graphic
, itd. W końcu funkcja 0x01 graphic
zeruje się w co najmniej jednym punkcie 0x01 graphic
należącym do najmniejszego przedziału zawierającego 0x01 graphic
. Wobec tego, że 0x01 graphic
, a 0x01 graphic
-sza pochodna wielomianu 0x01 graphic
wynosi 0x01 graphic
, mamy

0x01 graphic

Stąd

0x01 graphic

co kończy dowód.

Zwykle interesuje nas nie tyle błąd w ustalonym punkcie 0x01 graphic
, ale na całym przedziale 0x01 graphic
. Zakładając teraz, że przedział 0x01 graphic
jest domknięty, czyli

0x01 graphic

dla pewnych 0x01 graphic
, błąd ten będziemy mierzyć w normie jednostajnej (Czebyszewa). Dla funkcji ciągłej 0x01 graphic
, norma ta jest zdefiniowana jako

0x01 graphic

Niech 0x01 graphic
, gdzie 0x01 graphic
, będzie klasą funkcji

0x01 graphic

gdzie 0x01 graphic
. Mamy następujące twiedzenie.

Twierdzenie O najgorszym możliwym błędzie interpolacji w klasie

Załóżmy, że każdą funkcję 0x01 graphic
aproksymujemy jej wielomianem interpolacyjnym 0x01 graphic
opartym na 0x01 graphic
węzłach 0x01 graphic
. Wtedy maksymalny błąd takiej aproksymacji wynosi

0x01 graphic

Dowód. Oszacowanie górne wynika bezpośrednio z lematu o postaci błędu interpolacji, bowiem dla 0x01 graphic
mamy

0x01 graphic

Z drugiej strony zauważmy, że dla wielomianu 0x01 graphic
mamy 0x01 graphic
oraz

0x01 graphic

co kończy dowód.

Zjawisko Rungego i dobór węzłów interpolacji

Rozważmy zadanie interpolacji funkcji

0x01 graphic

w 0x01 graphic
równoodległych węzłach na przedziale 0x01 graphic
. Okazuje się, że dla dużych wartości 0x01 graphic
, wielomian interpolacyjny ma poważne kłopoty z aproksymacją tej funkcji przy krańcach przedziału:

Zjawisko Rungego: interpolacja w 0x01 graphic
węzłach równoodległych dla 0x01 graphic

Z kolei wielomian oparty na węzłach Czebyszewa znacznie lepiej przybliża tę funkcję.

Zjawisko Rungego: interpolacja w węzłach równoodległych, kontra interpolacja w węzłach Czebyszewa

Rzeczywiście, węzły Czebyszewa zagęszczają się w pobliżu krańców odcinka.

Zjawisko Rungego: interpolacja w węzłach Czebyszewa

Wiąże się to z zachowaniem się samych wielomianów bazowych: wielomiany oparte na węzłach równoodległych właśnie silnie oscylują w pobliżu krańców przedziału (jasne: nasz wielomian jest wysokiego stopnia, musi mieć dużo zer, a z drugiej strony, jako wielomian wysokiego stopnia, chce szybko uciec do nieskończoności, dlatego "wije się" jak może). Natomiast wielomiany bazowe oparte na węzłach Czebyszewa są najspokojniejsze: wiją się, ale z umiarem, bo zagęszczone przy krańcach węzły skutecznie je "duszą".

Zauważmy, że błąd aproksymacji 0x01 graphic
w istotny sposób zależy od wyboru węzłów 0x01 graphic
. Naturalne jest więc teraz następujące pytanie: w których punktach 0x01 graphic
przedziału 0x01 graphic
należy obliczać wartości funkcji, aby błąd był minimalny? Problem ten sprowadza się oczywiście do minimalizacji wielkości 0x01 graphic
względem węzłów 0x01 graphic
.

Twierdzenie O optymalnym doborze węzłów. Błąd aproksymacji w klasie funkcji 0x01 graphic
jest minimalny gdy węzły interpolacji są zadane jako węzły Czebyszewa na 0x01 graphic
, tzn.

0x01 graphic

Ponadto, dla optymalnych węzłów 0x01 graphic
mamy

0x01 graphic

Dowód tego twierdzenia opiera się na własnościach pewnego ważnego ciągu wielomianów, który teraz przedstawimy.

Wielomiany Czebyszewa

0x08 graphic
Ciąg 0x01 graphic
wielomianów Czebyszewa (pierwszego rodzaju) zdefiniowany jest indukcyjnie jako

0x01 graphic

Zauważmy, że 0x01 graphic
jest wielomianem stopnia dokładnie 0x01 graphic
o współczynniku przy 0x01 graphic
równym 0x01 graphic
(0x01 graphic
). Ponadto wielomian 0x01 graphic
można dla 0x01 graphic
przedstawić w postaci

0x01 graphic

Pafnutij Lwowicz Czebyszew

Rzeczywiście, łatwo sprawdzić, że jest to prawdą dla 0x01 graphic
. Stosując podstawienie 0x01 graphic
, 0x01 graphic
, oraz wzór na sumę cosinusów otrzymujemy dla 0x01 graphic
:

0x01 graphic

co jest równoważne formule rekurencyjnej dla 0x01 graphic
.

Kilka pierwszych wielomianów Czebyszewa na odcinku 0x01 graphic

Ze wzoru 0x01 graphic
wynikają również inne ważne własności wielomianów Czebyszewa. Norma wielomianu Czebyszewa na 0x01 graphic
wynosi

0x01 graphic

i jest osiągana w 0x01 graphic
punktach tego przedziału równych

0x01 graphic

przy czym 0x01 graphic
.

W końcu, 0x01 graphic
-ty wielomian Czebyszewa 0x01 graphic
ma dokładnie 0x01 graphic
pojedynczych zer w 0x01 graphic
równych

0x01 graphic

Miejsca zerowe wielomianu Czebyszewa będziemy nazywać węzłami Czebyszewa. Konsekwencją wymienionych własności jest następująca własność ekstremalna wielomianów Czebyszewa.

Przez 0x01 graphic
oznaczymy klasę wielomianów stopnia 0x01 graphic
o współczynniku wiodącym równym 0x01 graphic
, tzn.

0x01 graphic

Twierdzenie O minimaksie. Niech 0x01 graphic
. W klasie 0x01 graphic
minimalną normę jednostajną na przedziale 0x01 graphic
ma wielomian 0x01 graphic
, tzn.

0x01 graphic

Wielomian stopnia 9 oparty na węzłach Czebyszewa kontra oparty na węzłach równoodległych. Zwróć uwagę na wielkie oscylacje tego drugiego pry końcach odcinka.

Możemy teraz przeprowadzić dowód twierdzenia o optymalnym doborze węzłów. Dowód wynika teraz bezpośrednio z twierdzenia o minimaksie. Zauważmy bowiem, że wielomian 0x01 graphic
jest w klasie 0x01 graphic
. Stąd dla 0x01 graphic
optymalnymi węzłami są zera 0x01 graphic
wielomianu Czebyszewa, przy których

0x01 graphic

Jeśli przedział 0x01 graphic
jest inny niż 0x01 graphic
, należy dokonać liniowej zamiany zmiennych tak, aby przeszedł on na 0x01 graphic
. Bezpośrednie sprawdzenie pokazuje, że w klasie 0x01 graphic
minimalną normę Czebyszewa na przedziale 0x01 graphic
ma wielomian

0x01 graphic

Stąd

0x01 graphic

i węzły 0x01 graphic
są optymalne.

Wielomiany Czebyszewa znajdują bardzo wiele, czasem zaskakujących, zastosowań w różnych działach numeryki, m.in. w konstrukcji metod iteracyjnych rozwiązywania równań liniowych.

Równie interesujący jest fakt, że wielomian interpolacyjny oparty na węzłach Czebyszewa jest prawie optymalnym przybliżeniem wielomianowym zadanej funkcji:

Twierdzenie Jacksona, o prawie optymalnej interpolacji w węzłach Czebyszewa. Dla 0x01 graphic
, wielomian interpolacyjny 0x01 graphic
stopnia co najwyżej 0x01 graphic
, oparty na węzłach Czebyszewa, spełnia

0x01 graphic

gdzie 0x01 graphic
jest wielomianem stopnia co najwyżej 0x01 graphic
, najlepiej aproksymującym 0x01 graphic
w sensie normy jednostajnej.

Jeśli więc 0x01 graphic
, to wielomian oparty na węzłach Czebyszewa jest co najwyżej 3.02 razy, a gdy 0x01 graphic
--- maksymalnie 4 razy gorszy od optymalnego. Można więc powiedzieć, że jest prawie optymalny.

Literatura



Wyszukiwarka

Podobne podstrony:
MN 02 Row Nielin, metody numeryczne
Interpolacje na Metody Numeryczne w VBA excel
MN energetyka zadania od wykładowcy 09-05-14, STARE, Metody Numeryczne, Część wykładowa Sem IV
Cichy B Metody numeryczne, mn 09
Metody numeryczne PDF, MN calk 09
MN energetyka zadania od wykładowcy 09-05-14, STARE, Metody Numeryczne, Część wykładowa Sem IV
Metody numeryczne PDF, MN macierze 01 1
Metody numeryczne PDF, MN raphson 11
metody numeryczne - interpolacja, Nauka i Technika, Informatyka, Programowanie
02 Wybrane metody numeryczne (aproksymacja funkcji, rozwiazy
Cichy B Metody numeryczne, mn 06
Cichy B Metody numeryczne, mn 08
Metody numeryczne PDF, MN mnk1 06
Metody numeryczne PDF, MN inter 05
SCIAGA METODY NUMERYCZNE testy 1-8, Mechatronika, Semestr IV, Metody numeryczne, opracowanie MN, TES
MN 07 Uklady Row Lin 2, metody numeryczne
Metody numeryczne PDF, MN mnk2 07

więcej podobnych podstron