obliczenia_naukowe/l1/sprawozdanie.tex
2024-11-11 00:21:17 +01:00

536 lines
21 KiB
TeX
Raw Permalink Blame History

This file contains invisible Unicode characters

This file contains invisible Unicode characters that are indistinguishable to humans but may be processed differently by a computer. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

\documentclass[a4paper]{article}
\usepackage[T1]{fontenc}
\usepackage[polish]{babel}
\usepackage[utf8]{inputenc}
\usepackage{amsmath}
\usepackage{adjustbox}
\usepackage{longtable}
\usepackage{siunitx}
\title{Obliczenia naukowe Lista 1}
\author{Jacek Poziemski 272389}
\date{27.10.2024}
\begin{document}
\maketitle
\textbf{Zad. 1}
Wyznaczanie iteracyjnie epsilonów maszynowych, liczb maszynowych
$eta$ i maksymalnych wartości dla wszystkich dostępnych typów
zmiennopozycyjnych zgodnych ze standardem IEEE 754.
\vspace{0.5cm}
\textbf{Epsilon maszynowy}
\vspace{0.25cm}
Epsilon maszynowy to najmniejsza wartość dla danego typu
zmiennopozycyjnego, dla której spełniona jest następująca
nierówność: $1.0 + \varepsilon > 1.0$. Innymi słowami,
jest to wartość, która po dodaniu do liczby da nam następną
liczbę po niej. Epsilon maszynowy jest odwrotnie proporcjonalny
do precyzji typu zmiennopozycyjnego, w którym operujemy.
\vspace{0.25cm}
Epsilon liczę zaczynając od wartości $1.0$ dla danego typu i
dzieląc ją przez $2.0$ dopóki $1.0 + \frac{\varepsilon}{2.0} > 1.0$.
Podzielenie $\varepsilon$ przez $2.0$ w warunku pętli powoduje
wyjście z pętli po otrzymaniu najmniejszej wartości,
dla której warunek jest spełniony.
\vspace{0.25cm}
Porównanie wartości $\varepsilon$ liczonych iteracyjnie,
za pomocą funkcji $eps()$ z Julii oraz z headera
$\langle float.h \rangle$ biblioteki standardowej języka C:
\begin{table}[h]
\begin{adjustbox}{max width=\textwidth}
\begin{tabular}{|c|c|c|c|}
\hline
Typ zmiennopozycyjny & Wartość wyznaczona iteracyjnie & $eps()$ & $\langle float.h \rangle$ \\ \hline
Float16 & 0.000977 & 0.000977 &\\
Float32 & 1.1920929e-7 & 1.1920929e-7 & 1.1920929e-07 \\
Float64 & 2.220446049250313e-16 & 2.220446049250313e-16 & 2.220446049250313e-16 \\
\hline
\end{tabular}
\end{adjustbox}
\end{table}
Brak wartości Float16 dla headera $\langle float.h \rangle$
wynika z braku 16bitowego typu zmiennopozycyjnego w standardzie języka C.
\vspace{0.25cm}
Wartości wyliczone iteracyjnie są identyczne do tych z
bibliotek standardowych języków Julia i C, więc jeśli ufamy
twórcom tych języków, metoda iteracyjna jest wystarczająco dokładna.
\vspace{0.5cm}
\textbf{Liczba maszynowa eta}
\vspace{0.25cm}
Liczba maszynowa eta to najmniejsza zdenormalizowana
wartość dla danego typu zmiennopozycyjnego większa od $0.0$.
\vspace{0.25cm}
Liczbę eta liczę podobnie do liczby maszynowej epsilon
zaczynam od wartości $1.0$ i dzielę ją przez $2.0$ dopóki
warunek $\frac{\eta}{2.0} > 0.0$ jest spełniony.
Tak samo jak poprzednio podzielenie $\eta$ przez $2.0$
powoduje wyjście z pętli po otrzymaniu najmniejszej
wartości, dla której warunek jest spełniony.
Porównanie wartości $\eta$ liczonych iteracyjnie, za
pomocą $nextfloat(0.0)$ z Julii oraz wartościami zwróconymi
przez funkcję $floatmin()$ dla danego typu zmiennopozycyjnego:
\begin{table}[h]
\begin{adjustbox}{max width=\textwidth}
\begin{tabular}{|c|c|c|c|}
\hline
Typ zmiennopozycyjny & Wartość wyznaczona iteracyjnie & $nextfloat(0.0)$ & $floatmin()$ \\ \hline
Float16 & 6.0e-8 & 6.0e-8 & 6.104e-5 \\
Float32 & 1.0e-45 & 1.0e-45 & 1.1754944e-38 \\
Float64 & 5.0e-324 & 5.0e-324 & 2.2250738585072014e-308 \\
\hline
\end{tabular}
\end{adjustbox}
\end{table}
Wartości wyliczone iteracyjnie są identyczne do tych
zwróconych przez funkcję $nextfloat(0.0)$, ale $floatmin()$
zwraca wartości o wiele większe. Jest tak dlatego, że ta
funkcja zwraca wartości znormalizowane, a metoda iteracyjna
i funkcja $nextfloat(0.0)$ zwracają wartości zdenormalizowane
takie, których cecha składa się z samych bitów 0.
\vspace{0.5cm}
\textbf{Maksymalne wartości}
\vspace{0.5cm}
Wartość maksymalną dla danego typu zmiennopozycyjnego
liczę zaczynając od $1.0$ dla tego typu. Następnie
mnożę tą wartość przez $2.0$ dopóki $max \cdot 2.0 < \infty$
(sprawdzam tę nierówność korzystając z funkcji $isfinite()$).
Osiągnięta w ten sposób wartość niekoniecznie jest maksymalna,
żeby to zmienić tworzę pomocniczą zmienną $x$, której na początku
przypisuję wartość $\frac{max}{2.0}$. Następnie, dopóki
$max + x < \infty$ dodaję $x$ do $max$ i dzielę $x$ przez $2.0$.
\vspace{0.25cm}
W przypadku gdyby w ostatniej iteracji pierwszej pętli $max$
nie był maksymalny, ale $max \cdot 2.0$ był nieskończony,
zostało nam co najwyżej $max \cdot 2.0 - \varepsilon$ do
"`dopełnienia"'. Zaczynamy to "`dopełnianie"' od $\frac{max}{2.0}$
wiemy, że $max + max = max \cdot 2.0$ będzie nieskończone.
W każdej iteracji dodajemy $x$ do $max$ i dzielimy $x$ przez $2.0$.
W ten sposób "`dopełnimy"' $max$ do końca i otrzymamy
rzeczywistą wartość maksymalną.
\vspace{0.25cm}
Porównanie wartości maksymalnych liczonych iteracyjnie,
za pomocą funkcji $floatmax()$ z Julii oraz z headera
$\langle float.h \rangle$ biblioteki standardowej języka C:
\begin{table}[h]
\begin{adjustbox}{max width=\textwidth}
\begin{tabular}{|c|c|c|c|}
\hline
Typ zmiennopozycyjny & Wartość wyznaczona iteracyjnie & $floatmax()$ & $\langle float.h \rangle$ \\ \hline
Float16 & 6.55e4 & 6.55e4 &\\
Float32 & 3.4028235e38 & 3.4028235e38 & 3.4028235e+38 \\
Float64 & 1.7976931348623157e308 & 1.7976931348623157e308 & 1.7976931348623157e+308 \\
\hline
\end{tabular}
\end{adjustbox}
\end{table}
Brak wartości Float16 dla headera $\langle float.h \rangle$
wynika z braku 16bitowego typu zmiennopozycyjnego w standardzie języka C.
\vspace{0.25cm}
Wartości wyznaczone iteracyjnie są identyczne do tych
zwróconych przez $floatmax()$ jak i tych z headera
$\langle float.h \rangle$ z biblioteki standardowej języka C.
\vspace{0.5cm}
\textbf{Zad. 2}
Eksperymentalne sprawdzenie słuszności twierdzenia Kahana dla
wszystkich typów zmiennopozycyjnych dostępnych w języku Julia.
\vspace{0.25cm}
Twierdzenie Kahana mówi, że epsilon maszynowy można
uzyskać w następujący sposób: $\varepsilon = 3(4/3 - 1) - 1$.
Aby sprawdzić wyniki korzystam z funkcji
$eps()$ z biblioteki standardowej Julii.
\vspace{0.25cm}
Porównanie wartości epsilonu maszynowego liczonego za
pomocą twierdzenia Kahana oraz za pomocą funkcji $eps()$:
\begin{table}[h]
\begin{adjustbox}{max width=\textwidth}
\begin{tabular}{|c|c|c|}
\hline
Typ zmiennopozycyjny & Twierdzenie Kahana & $eps()$ \\ \hline
Float16 & -0.000977 & 0.000977 \\
Float32 & 1.1920929e-7 & 1.1920929e-7 \\
Float64 & -2.220446049250313e-16 & 2.220446049250313e-16 \\
\hline
\end{tabular}
\end{adjustbox}
\end{table}
Dla typu $Float32$ wartości są identyczne, dla typów
$Float16$ i $Float64$ wartości różnią się znakiem
twierdzenie Kahana daje nam wartości ujemne. Wynika
to z definicji typów zmiennoprzecinkowych i rozwinięcia
liczby $4/3$ binarnie $1.(10)$. Typy $Float16$, $Float32$ i
$Float64$ mają odpowiednio 10, 23 i 52 bity mantysy.
Rozwijając $4/3$, na indeksach parzystych dostajemy bit 0,
na nieparzystych bit 1. W typach o parzystej ilości bitów
w mantysie ($Float16$ i $Float64$) na ostatnim indeksie będzie 0.
Mnożąc wartość $4/3 - 1$ przez $3$ otrzymamy wartość minimalnie
mniejszą od $1$ Przez co wynik jest ujemny. Z kolei dla typu
$Float32$ ostatnim bitem będzie 1, co da nam wartość dodatnią.
\vspace{0.25cm}
Twierdzenie Kahana poprawnie liczy epsilon maszynowy,
jednak trzeba pamiętać o wzięciu modułu z wyniku.
\vspace{0.5cm}
\textbf{Zad. 3}
Eksperymentalne sprawdzenie rozmieszczenia liczb
zmiennopozycyjnych w arytmetyce $Float64$ w języku Julia.
\vspace{0.25cm}
Aby sprawdzić, czy liczby w przedziale $[1, 2]$
równomiernie rozmieszczone można zrobić pętlę,
w której od $1.0$ do $2.0$ jedną liczbę iterujemy
dodając do niej krok $\delta = 2^{-52}$, drugą korzystając
z funkcji $nextfloat()$ jednak jest to bardzo nieefektywna metoda.
\vspace{0.25cm}
Lepszą metodą jest sprawdzenie pierwszej i ostatniej
liczby z przedziału, w naszym przypadku $nextfloat(1.0)$
i $nextfloat(2.0)$. Sprawdzamy cechę obu jeśli są różne,
wyklucza to równomierne rozmieszczenie liczb. W przeciwnym
przypadku możemy obliczyć $\delta$ w następujący sposób:
\[2^{cecha - 1023} \cdot 2^{-52}\]
gdzie 1023 to przesunięcie cechy typu $Float64$, a 52 to rozmiar mantysy.
\begin{table}[h]
\begin{adjustbox}{max width=\textwidth}
$\begin{array}{|c|c|}
\hline
\text{Przedział} & \delta \\ \hline
\text{[1.0, 2.0]} & \num{2.220446049250313e-16} = 2^{-52} \\
\text{[0.5, 1.0]} & \num{1.1102230246251565e-16} = 2^{-53} \\
\text{[2.0, 4.0]} & \num{4.440892098500626e-16} = 2^{-51} \\
\hline
\end{array}$
\end{adjustbox}
\end{table}
Dla przedziału [1.0, 2.0] rzeczywiście $\delta = 2^{-52}$,
z kolei $\delta$ dla przedziałów [0.5, 1.0] i [2.0, 4.0]
pokazuje, że jest ona zależna od cechy liczb z zakresu
nad którym się skupiamy. Im bliżej zera jesteśmy tym mniejsza
jest $\delta$ zwiększa się dokładność liczb. Analogicznie
oddalając się od zera tracimy precyzję części ułamkowej
przy coraz większych liczbach prawdopodobnie bardziej obchodzi nas część całkowita.
\vspace{0.5cm}
\textbf{Zad. 4}
Eksperymentalne znalezienie liczby w arytmetyce $Float64$ w
przedziale $(1, 2)$ takiej, że $x \cdot (1/x) \neq 1$ oraz
znalezienie najmniejszej takiej liczby.
\vspace{0.25cm}
Implementacja zadania opiera się na funkcji $nextfloat()$ z
języka Julia. Zaczynamy od $nextfloat(1.0)$, iterujemy do
$2.0$. W każdej iteracji liczymy wynik $x \cdot (1/x)$,
jeśli jest różny od $1.0$ zwracamy $x$ jako wynik.
Oczywiście jest to również najmniejsza taka liczba.
\vspace{0.25cm}
Najmniejsza taka liczba to $1.000000057228997$. Błędny wynik
wyrażenia, które zawsze powinno zwracać $1.0$ jest spowodowany
oczywiście skończoną mantysą w liczbach zmiennoprzecinkowych
przez co tracimy precyzję ale również kolejnością wykonywania
działań. Gdyby zamienić $x \cdot (1/x)$ na $(x \cdot 1) / x$
zawsze dostalibyśmy odpowiedni wynik $x \cdot 1$ nie
zmieniłoby wartości liczby, a $x / x$ zakładając dobrze
napisany interpreter lub kompilator zawsze zwróci $1.0$
\vspace{0.5cm}
\textbf{Zad. 5}
Eksperymentalne obliczanie iloczynu skalarnego dwóch wektorów:
\begin{gather*}
x = [2.718281828, -3.141592654, 1.414213562, 0.5772156649, 0.3010299957] \\
y = [1486.2497, 878366.9879, -22.37492, 4773714.647, 0.000185049]
\end{gather*}
na 4 sposoby:
\vspace{0.25cm}
(a) w przód od pierwszego do ostatniego elementu
\vspace{0.25cm}
(b) w tył od ostatniego do pierwszego elementu
\vspace{0.25cm}
(c) od największego do najmniejszego liczenie częściowych
sum elementów dodatnich i ujemnych w porządku malejącym względem
wartości absolutnej liczb, następnie dodanie do siebie sum częściowych
\vspace{0.25cm}
(d) od największego do najmniejszego liczenie częściowych
sum elementów dodatnich i ujemnych w porządku rosnącym względem
wartości absolutnej liczb, następnie dodanie do siebie sum częściowych
\vspace{0.25cm}
używając typów zmiennopozycyjnych $Float32$ i $Float64$.
\vspace{0.25cm}
Oczekiwanym wynikiem iloczynu skalarnego na tych
dwóch wektorach jest $-1.00657107000000 \cdot 10^{-11}$.
\vspace{0.25cm}
Porównanie wyników iloczynu skalarnego przeprowadzonego
za pomocą powyższych czterech sposobów:
\begin{table}[h]
\begin{adjustbox}{max width=\textwidth}
\begin{tabular}{|c|c|c|}
\hline
Sposób & Wynik dla $Float32$ & Wynik dla $Float64$ \\ \hline
(a) & -0.4999442994594574 & 1.0251881368296672e-10 \\
(b) & -0.454345703125 & -1.5643308870494366e-10 \\
(c) & -0.5 & 0.0 \\
(d) & -0.5 & 0.0 \\
\hline
\end{tabular}
\end{adjustbox}
\end{table}
Żaden ze sposobów nie dał wyniku bliskiego oczekiwanemu,
choć sposób (b) w arytmetyce $Float64$ dał wynik dość
bliski prawdzie. Patrząc na wyniki korzystające z różnych
sposobów widzimy również, że kolejność wykonywania działań
może znacznie wpłynąć na poprawność wyniku.
\vspace{0.5cm}
\textbf{Zad. 6}
Policzenie w języku Julia w arytmetyce $Float64$ wartości funkcji:
\[f(x) = \sqrt{x^2 + 1} - 1\]
\[g(x) = x^2 / (\sqrt{x^2 + 1} + 1)\]
dla argumentów $x = 8^{-1}, 8^{-2}, 8^{-3}, ...$.
\vspace{0.25cm}
Porównanie wyników dla obu funkcji:
\newpage
\begin{table}[h]
\begin{adjustbox}{max width=\textwidth}
\begin{tabular}{|c|c|c|}
\hline
$x$ & $f(x)$ & $g(x)$ \\ \hline
$8^{-1}$ & 0.0077822185373186414 & 0.0077822185373187065 \\
$8^{-2}$ & 0.00012206286282867573 & 0.00012206286282875901 \\
$8^{-3}$ & 1.9073468138230965e-6 & 1.907346813826566e-6 \\
$8^{-4}$ & 2.9802321943606103e-8 & 2.9802321943606116e-8 \\
$8^{-5}$ & 4.656612873077393e-10 & 4.6566128719931904e-10 \\
$8^{-6}$ & 7.275957614183426e-12 & 7.275957614156956e-12 \\
$8^{-7}$ & 1.1368683772161603e-13 & 1.1368683772160957e-13 \\
$8^{-8}$ & 1.7763568394002505e-15 & 1.7763568394002489e-15 \\
$8^{-9}$ & 0.0 & 2.7755575615628914e-17 \\
$8^{-10}$ & 0.0 & 4.336808689942018e-19 \\
$8^{-11}$ & 0.0 & 6.776263578034403e-21 \\
$8^{-12}$ & 0.0 & 1.0587911840678754e-22 \\
$8^{-13}$ & 0.0 & 1.6543612251060553e-24 \\
$8^{-14}$ & 0.0 & 2.5849394142282115e-26 \\
$8^{-15}$ & 0.0 & 4.0389678347315804e-28 \\
$8^{-16}$ & 0.0 & 6.310887241768095e-30 \\
$8^{-17}$ & 0.0 & 9.860761315262648e-32 \\
$8^{-18}$ & 0.0 & 1.5407439555097887e-33 \\
$8^{-19}$ & 0.0 & 2.407412430484045e-35 \\
$8^{-20}$ & 0.0 & 3.76158192263132e-37 \\
\hline
\end{tabular}
\end{adjustbox}
\end{table}
Już przy $x = 8^{-9}$, $f(x)$ przyjmuje wartość $0.0$,
którą zachowuje przy mniejszych wartościach. Przy
wartościach większych od $8^{-9}$ obie funkcje mają
podobne wartości. $f(x)$ szybciej traci precyzję i
spada do $0.0$ przez odejmowanie $1.0$ od pierwiastka
im mniejsza wartość $x$ tym bliżej $x^2 + 1$ jest $1.0$.
Przez pierwiastek całe wyrażenie zbliża się jeszcze bardziej
do $1.0$, od którego również odejmujemy $1.0$.
Przekształcenie $f(x)$ w $g(x)$ zwalcza ten problem
poprzez pozbycie się problematycznego odejmowania,
dzięki czemu tracimy mniej precyzji wykonując działania.
\vspace{0.5cm}
\textbf{Zad. 7}
Obliczanie przybliżonej wartości pochodnej funkcji
$f(x) = \sin{x} + \cos{3x}$ w punkcie $x_0 = 1$ oraz
błędów $|f^{'}(x_0) - \widetilde{f^{'}}(x_0)|$ korzystając ze wzoru
\[f^{'}(x_0) \approx \widetilde{f^{'}}(x_0) = \frac{f(x_0 + h) - f(x_0)}{h}\]
dla $h = 2^{-n}$, gdzie $n = 0, 1, 2, ..., 54$.
\vspace{0.25cm}
Dokładna pochodna $f(x)$ to $f^{'}(x) = \cos{x} - 3\sin{3x}$
\vspace{0.25cm}
Dokładna wartość pochodnej dla $x_0 = 1$ to $0.11694228168853815$.
\begin{longtable}{|c|c|c|c|c|}
\hline
$n$ & $h$ & $1 + h$ & $f^{'}(x_0)$ & błąd \\ \hline
$0$ & 1.0 & 2.0 & 2.0179892252685967 & 1.9010469435800585 \\
$1$ & 0.5 & 1.5 & 1.8704413979316472 & 1.753499116243109 \\
$2$ & 0.25 & 1.25 & 1.1077870952342974 & 0.9908448135457593 \\
$3$ & 0.125 & 1.125 & 0.6232412792975817 & 0.5062989976090435 \\
$4$ & 0.0625 & 1.0625 & 0.3704000662035192 & 0.253457784514981 \\
$5$ & 0.03125 & 1.03125 & 0.24344307439754687 & 0.1265007927090087 \\
$6$ & 0.015625 & 1.015625 & 0.18009756330732785 & 0.0631552816187897 \\
$7$ & 0.0078125 & 1.0078125 & 0.1484913953710958 & 0.03154911368255764 \\
$8$ & 0.00390625 & 1.00390625 & 0.1327091142805159 & 0.015766832591977753 \\
$9$ & 0.001953125 & 1.001953125 & 0.1248236929407085 & 0.007881411252170345 \\
$10$ & 0.0009765625 & 1.0009765625 & 0.12088247681106168 & 0.0039401951225235265 \\
$11$ & 0.00048828125 & 1.00048828125 & 0.11891225046883847 & 0.001969968780300313 \\
$12$ & 0.000244140625 & 1.000244140625 & 0.11792723373901026 & 0.0009849520504721099 \\
$13$ & 0.0001220703125 & 1.0001220703125 & 0.11743474961076572 & 0.0004924679222275685 \\
$14$ & 6.103515625e-5 & 1.00006103515625 & 0.11718851362093119 & 0.0002462319323930373 \\
$15$ & 3.0517578125e-5 & 1.000030517578125 & 0.11706539714577957 & 0.00012311545724141837 \\
$16$ & 1.52587890625e-5 & 1.0000152587890625 & 0.11700383928837255 & 6.155759983439424e-5 \\
$17$ & 7.62939453125e-6 & 1.0000076293945312 & 0.11697306045971345 & 3.077877117529937e-5 \\
$18$ & 3.814697265625e-6 & 1.0000038146972656 & 0.11695767106721178 & 1.5389378673624776e-5 \\
$19$ & 1.9073486328125e-6 & 1.0000019073486328 & 0.11694997636368498 & 7.694675146829866e-6 \\
$20$ & 9.5367431640625e-7 & 1.0000009536743164 & 0.11694612901192158 & 3.8473233834324105e-6 \\
$21$ & 4.76837158203125e-7 & 1.0000004768371582 & 0.1169442052487284 & 1.9235601902423127e-6 \\
$22$ & 2.384185791015625e-7 & 1.000000238418579 & 0.11694324295967817 & 9.612711400208696e-7 \\
$23$ & 1.1920928955078125e-7 & 1.0000001192092896 & 0.11694276239722967 & 4.807086915192826e-7 \\
$24$ & 5.960464477539063e-8 & 1.0000000596046448 & 0.11694252118468285 & 2.394961446938737e-7 \\
$25$ & 2.9802322387695312e-8 & 1.0000000298023224 & 0.116942398250103 & 1.1656156484463054e-7 \\
$26$ & 1.4901161193847656e-8 & 1.0000000149011612 & 0.11694233864545822 & 5.6956920069239914e-8 \\
$27$ & 7.450580596923828e-9 & 1.0000000074505806 & 0.11694231629371643 & 3.460517827846843e-8 \\
$28$ & 3.725290298461914e-9 & 1.0000000037252903 & 0.11694228649139404 & 4.802855890773117e-9 \\
$29$ & 1.862645149230957e-9 & 1.0000000018626451 & 0.11694222688674927 & 5.480178888461751e-8 \\
$30$ & 9.313225746154785e-10 & 1.0000000009313226 & 0.11694216728210449 & 1.1440643366000813e-7 \\
$31$ & 4.656612873077393e-10 & 1.0000000004656613 & 0.11694216728210449 & 1.1440643366000813e-7 \\
$32$ & 2.3283064365386963e-10 & 1.0000000002328306 & 0.11694192886352539 & 3.5282501276157063e-7 \\
$33$ & 1.1641532182693481e-10 & 1.0000000001164153 & 0.11694145202636719 & 8.296621709646956e-7 \\
$34$ & 5.820766091346741e-11 & 1.0000000000582077 & 0.11694145202636719 & 8.296621709646956e-7 \\
$35$ & 2.9103830456733704e-11 & 1.0000000000291038 & 0.11693954467773438 & 2.7370108037771956e-6 \\
$36$ & 1.4551915228366852e-11 & 1.000000000014552 & 0.116943359375 & 1.0776864618478044e-6 \\
$37$ & 7.275957614183426e-12 & 1.000000000007276 & 0.1169281005859375 & 1.4181102600652196e-5 \\
$38$ & 3.637978807091713e-12 & 1.000000000003638 & 0.116943359375 & 1.0776864618478044e-6 \\
$39$ & 1.8189894035458565e-12 & 1.000000000001819 & 0.11688232421875 & 5.9957469788152196e-5 \\
$40$ & 9.094947017729282e-13 & 1.0000000000009095 & 0.1168212890625 & 0.0001209926260381522 \\
$41$ & 4.547473508864641e-13 & 1.0000000000004547 & 0.116943359375 & 1.0776864618478044e-6 \\
$42$ & 2.2737367544323206e-13 & 1.0000000000002274 & 0.11669921875 & 0.0002430629385381522 \\
$43$ & 1.1368683772161603e-13 & 1.0000000000001137 & 0.1162109375 & 0.0007313441885381522 \\
$44$ & 5.684341886080802e-14 & 1.0000000000000568 & 0.1171875 & 0.0002452183114618478 \\
$45$ & 2.842170943040401e-14 & 1.0000000000000284 & 0.11328125 & 0.003661031688538152 \\
$46$ & 1.4210854715202004e-14 & 1.0000000000000142 & 0.109375 & 0.007567281688538152 \\
$47$ & 7.105427357601002e-15 & 1.000000000000007 & 0.109375 & 0.007567281688538152 \\
$48$ & 3.552713678800501e-15 & 1.0000000000000036 & 0.09375 & 0.023192281688538152 \\
$49$ & 1.7763568394002505e-15 & 1.0000000000000018 & 0.125 & 0.008057718311461848 \\
$50$ & 8.881784197001252e-16 & 1.0000000000000009 & 0.0 & 0.11694228168853815 \\
$51$ & 4.440892098500626e-16 & 1.0000000000000004 & 0.0 & 0.11694228168853815 \\
$52$ & 2.220446049250313e-16 & 1.0000000000000002 & -0.5 & 0.6169422816885382 \\
$53$ & 1.1102230246251565e-16 & 1.0 & 0.0 & 0.11694228168853815 \\
$54$ & 5.551115123125783e-17 & 1.0 & 0.0 & 0.11694228168853815 \\
\hline
\end{longtable}
Najdokładniejszy wynik jest dla $n = 28$, $h = 2^{-28}$,
dokładnie po środku wartości testowanych (jeśli zignorujemy
54, dla którego wynik jest identyczny do 53). Idąc w górę lub w dół
błąd jest coraz większy, dla $n = 53$ dochodząc do wartości pochodnej.
\vspace{0.25cm}
Wzięcie zbyt dużego $h$ (zbyt małego $n$) powoduje ogromny
błąd z definicji pochodnej, którą jest granica wzoru
wykorzystanego w tym zadaniu przy h dążącym do 0. Przy coraz
mniejszych $h$ (coraz większych $n$) na papierze zbliżamy się
do poprawnej wartości, ale w rzeczywistości ograniczona
precyzja liczb zmiennopozycyjnych w standardzie
IEEE 754 powoduje zwiększanie się błędu, lecz nie do
takich rozmiarów jak dla pierwszych kilku $n$.
Maksymalny błąd po przekroczeniu $n = 28$ otrzymujemy
przy $n = 53$, gdzie $h < \varepsilon$ w precyzji $Float64$. Wtedy
\[x_0 + h = x_0 = 1\]
a we wzorze na pochodną
\[\frac{f(x_0 + h) - f(x_0)}{h} = \frac{f(x_0) - f(x_0)}{h} = \frac{0}{h} = 0\]
zatem dla $n \geq 53$ błąd będzie zawsze równy dokładnej wartości pochodnej.
\vspace{0.5cm}
\textbf{Wnioski z listy}
Liczby zmiennopozycyjne w standardzie IEEE 754 mają skończoną
precyzję, co przy liczbach bliskich zeru powoduje błędy. Pisząc
kod operujący na tych liczbach musimy pamiętać o tym fakcie
i wybrać odpowiedni typ do trzymania danych lub zastanowić
się nad trzymaniem ich w typach całkowitych, np. waluty zamiast
w liczbie zmiennopozycyjnej trzymać w dwóch liczbach całkowitych.
\end{document}