VI. Elementy Metod Numerycznych
6.1 Podstawowe pojęcia związane z obliczeniami numerycznymi
6.1.1 Numeryczna reprezentacja liczb
Liczby w komputerze są reprezentowane przez skończoną liczbę cyfr. Sposób reprezentacji liczby zależy od jej typu np. liczby całkowite (typu integer) są przedstawiane w sposób stałopozycyjny, a rzeczywiste (typu real) w sposób zmiennopozycyjny.
Dowolną liczbę całkowitą k możemy przedstawić w postaci rozwinięcia dwójkowego
(en≠0 dla k≠0)
gdzie s=±1, a ei=0 lub 1 dla i=0,2,..,n. W komputerze na reprezentację liczby przeznacza się skończoną liczbę bitów np. d+1 (jeden bit na znak). Jeśli n<d to mówimy, że liczba k jest reprezentowana w rozpatrywanej arytmetyce. W ten sposób mogą być reprezentowane liczby całkowite z przedziału [-2d+1;2d-1]. Zakładając, że argumenty i wynik dodawania, odejmowania i mnożenia liczb całkowitych są reprezentowane to powyższe działania są wykonywane dokładnie.
Dowolną liczbę rzeczywistą x≠0 możemy przedstawić w postaci
x=s*2cm
gdzie s=±1, c - liczba rzeczywista zwana cechą, m - liczba rzeczywista zwana mantysą (m∈<0,5;1)). Dla x≠0 przedstawienie jest jednoznaczne. Cechę zapisujemy w sposób stałopozycyjny na d-t bitach. Pozostałe t bitów przeznaczone jest na mantysę. tak więc zapamiętujemy tylko t pierwszych bitów z na ogół nieskończonego rozwinięcia mantysy. (Jest to przyczyna niedokładności obliczeń na liczbach rzeczywistych). Jeżeli mt jest zaokrągloną do t miejsc po przecinku mantysą to łatwo widać, że
Zero jest zazwyczaj reprezentowane przez słowo o wszystkich bitach równych zero. Reprezentację zmiennopozycyjną liczby x będziemy oznaczali przez rd(x).
rd(x)=s*2cmt
Zakładając, że x≠0 prawdziwy jest poniższy wzór
albo inaczej rd(x)=x(1+ε).
Realizację działań w arytmetyce zmiennopozycyjnej będziemy oznaczali skrótem fl(ang. floating point)
6.1.2 Uwarunkowanie zadania
Zajmijmy się problemem polegającym na obliczeniu wyniku w=(w1,w2,..,wm) dla danych d(d1,d2,..,dn)
w=ϕ(d)
gdzie ϕ jest ciągłym odwzorowaniem ϕ:D0 ⊂ Rd → Rw, Rw i Rd są skończenie wymiarowymi przestrzeniami kartezjańskimi. Należy pamiętać, że do dyspozycji mamy tylko reprezentacjami rd(d1), rd(d2),..,rd(dn) danych początkowych.
Jeśli niewielkie względne zmiany danych początkowych zadania powodują duże względne zmiany jego rozwiązania, to zadanie takie nazywamy źle uwarunkowanym. Wielkości charakteryzujące wpływ zaburzeń danych początkowych na zaburzenia rozwiązania nazywamy wskaźnikami uwarunkowania zadania.
Przykład 6.1
Sprawdzić, czy problem obliczania funkcji f(x)=sin(x) jest dobrze uwarunkowany tzn. dla jakich argumentów x zadanie obliczenia sin(x) jest dobrze uwarunkowane.
Rozpatrzmy nieco zaburzony wynik dla nieco zaburzonych danych. f(x(1+ε1))(1+ε2)=f(x+xε1)(1+ε2). Korzystając ze wzoru Taylor'a możemy zauważyć, że f(x+Δx)=f(x)+f'(ξ)xε, gdzie ξ∈[x,x+xε].
, dlla f(x)≠0 i |ε|≤δ<<1.
zauważmy jeszcze, że |sin(x+Δx)-sin(x)|=|cos(x)Δx| + θ(Δx2)
czyli w naszym przypadku
dla x→0 błąd może być dowolnie duży, gdyż sin(x) może być dowolnie mały oraz dla x=kΠ zadanie jest źle uwarunkowane.
Proszę zapamiętać, że uwarunkowanie zadania jest cechą zadania niezależną od metody jego rozwiązywania.
6.1.3 Numeryczna poprawność algorytmów
Algorytm nazwiemy numerycznie poprawny jeśli obliczone rozwiązanie jest nieco zaburzonym rozwiązaniem( ale dokładnym) zadania o nieco zaburzonych danych. ”Nieco zaburzone” oznacza, że zaburzenie jest na poziomie błędu reprezentacji.
Innymi słowy algorytm A nazywamy numerycznie poprawnym w klasie zadań {ϕ,D} jeśli istnieją stałe Rd i Rw takie, że dla każdego d∈D istnieje d0 ∈ D0 takie, że
|d-d0|≤δdKd|d|
|ϕ(d0)-fl(A(d))|≤δwKw|ϕ(d0)|
Przykład 6.2
Zbadać, czy algorytm liczący iloczyn skalarny jest numerycznie poprawny.
Niech i oraz ; ; fl - skończona, zmiennopozycyjna arytmetyka to
Otrzymaliśmy dokładny wynik dla nieco zaburzonych danych, czyli algorytm jest numerycznie poprawny.
6.1.4 Algorytmy numerycznie stabilne
Algorytm A nazwiemy numerycznie stabilnym w klasie zadań {ϕ,D} jeśli istnieje stała K taka, że dla każdego d∈D zachodzi nierówność.
|ϕ(d)-fl(A(d))|≤K*P(d,ϕ)
gdzie P(d,ϕ) - optymalny poziom błędu rozwiązania ϕ(d) w arytmetyce fl.
Można sobie wyobrażać, że algorytm numerycznie stabilny to taki, który gwarantuje uzyskanie rozwiązania z błędem tego samego rządu, co optymalny poziom błędu rozwiązania danego zadania.
Łatwo można zauważyć, że każdy algorytm numerycznie poprawny jest algorytmem numerycznie stabilnym. Własności te nie są równoważne.
Zadania:
6.1.1
Zbadać, czy algorytm liczący iloczyn macierz przez wektor jest numerycznie poprawny.
6.1.2
Sprawdzić, jak błąd danych początkowych wpływa na wynik końcowy.
a) g(a,b)=a2+b2 b)g(a,b)=a2-b2
6.1.3
Zbadaj numeryczną poprawność algorytmu liczącego iloczyn skalarny dwóch wektorów n - elementowych.
6.2 Interpolacja
6.2.1 Interpolacja wielomianowa
Interpolacja Lagrange'a
Zadanie interpolacyjne Lagrange'a polega na znalezieniu dla danej funkcji f wielomianu Ln stopnia nie wyższego niż n, którego wartości w n+1 punktach xi są takie same jak wartości interpolowanej funkcji tzn.
Ln(xi)=f(xi), dla i=0,1,2,..,n
gdzie xi≠xj dla i≠j.
Wielomian Ln nazywamy wielomianem interpolacyjnym Lagrange'a funkcji f opartym na węzłach x0,x1,..,xn (parami różnych) jeśli
Wielomian interpolacyjny Ln(x) można również zapisać w postaci Newtona tzn.
,
gdzie pk(x) to wielomian zdefiniowany następująco
p0(x)=1
pk(x)=(x-x0)(x-x1)..(x-xk-1) dla k=1,2,..,n,
a współczynniki bk mają postać
i nazywane są ilorazami różnicowymi funkcji f opartymi na węzłach x0,x1,..,xk.
Iloraz różnicowy rzędu k funkcji f opartym na parami różnych węzłach x0,x1,..,xk (w których f jest określona) nazywamy wyrażenie f[x0,x1,..,xk] postaci
Iloraz różnicowy rzędu zerowego oparty na węźle x0 rozumiemy wartość funkcji w tym punkcie tzn. f[x0]=f(x0).
Tw.
Dla dowolnego układu parami różnych punktów x0,x1,..,xk należących do dziedziny funkcji f, zachodzi zależność
Przykład 6.3
Jak obliczyć pierwsze trzy współczynniki bk wielomianu interpolacyjnego opartego na węzłach 0, 1, 2, wiedząc, że f(0)=1, f(1)=2, f(2)=4?
Rozwiązanie:
Niech x0=0, x1=1, x2=2 to
Obliczymy je rekurencyjnie tzn.
x0 f(x0)=1
x1 f(x1)=2
x2 f(x2)=4
Tak więc nasz wielomian ma postać
L3(x)=1+1(x-0)+(x-0)(x-1)
Interpolacja Hermite'a
Danych jest k+1 różnych węzłów x0,x1,..,xk oraz liczby m0,m1,..,mk takie, że . Zadanie interpolacyjne Hermite'a polega na znalezieniu dla danej funkcji f wielomianu Hn stopnia nie wyższego niż n, spełniającego warunki
Hn(j)=f(j)(xj), dla i=0,1,..,k j=0,1,..,mi-1
Iloraz różnicowy oparty na i takich samych węzłach xl określamy następująco
Posiadając powyższą definicję możemy zdefiniować wielomian Hermite'a
,
gdzie bl=f[x0,m0;x1,m1;..;xi-1,mi-1;xi,j+1] bl=m0+m1+..+mi-1+j+1, a pl(x)=.
Jak można nietrudno współczynniki wielomianu interpolacyjnego Hermite'a oblicza się podobnie do współczynników wielomianu Lagrange'a.
zadania
6.2.1
Obliczyć wielomian Lagrange'a najmniejszego stopnia oparty na węzłach -1, 0, 2 wiedząc, że
a) f(-1)=0, f(0)=2 f(2)=1
b) f(-1)=-1, f(0)=0, f(2)=2
6.2.2
Mamy dany wielomian f(x)=x2+bx+c. Zbudować dla tego wielomianu ilorazy różnicowe: f[x,x0], f[x,x0,x1], f[x,x0,x1,x2]. Wykaż, że ostatni z ilorazów różnicowych jest tożsamościowo równy zero.
6.2.3
Obliczyć wartość ilorazów różnicowych funkcji f(x), która dla x=-1, 0, 2, 3, 5 przyjmuje wartości f(x)=1, 3, -11, -27, 103.
6.2.4
Zbudować wielomian interpolacyjny Hermite'a, który w punktach x=-1, 0, 2, 3, 5 przyjmuje wartości w(x)=1, 3, -11, -27, 103.
6.2.5
Zbudować wielomian interpolacyjny Hermite'a dla funkcji f(x) wiedząc, że f(0)=1, f'(0)=2, f(2)=3, f'(2)=3, f”(2)=4.
172