-
61. Data: 2014-08-05 08:00:26
Temat: Re: szybki logarytm
Od: Tomasz Kaczanowski <kaczus@dowyciecia_poczta.onet.pl>
W dniu 2014-08-05 01:15, slawek pisze:
> On Mon, 04 Aug 2014 20:58:28 +0200, Borneq <b...@a...hidden.pl>
> wrote:
>> Mnożenie jest bardzo szybkie, nie trzeba na intach, można na
> doublach, i
>> tak będzie szybko
>
> Super. Będzie szybko. I coraz mniej dokładnie.
Mnożenie jest operacją szybką i w zasadzie jedyna możliwość utraty
dokładności to konwersja inta na doubla (i na odwrót) - więc to zależy
czego oczekujemy.
--
Kaczus
http://kaczus.ppa.pl
-
62. Data: 2014-08-05 08:38:32
Temat: Re: szybki logarytm
Od: slawek <f...@f...com>
On Tue, 05 Aug 2014 08:00:26 +0200, Tomasz Kaczanowski
<kaczus@dowyciecia_poczta.onet.pl> wrote:
> Mnożenie jest operacją szybką i w zasadzie jedyna możliwość utraty
> dokładności to konwersja inta na doubla (i na odwrót) - więc to
zależy
Oblicz 1.000000001 do potęgi 100000 przez kolejne mnożenia i wzorem.
Oczywiście przez mnożeniem będzie szybciej jeżeli masz stabilicować
wartości dla wykładników od 1 do 100000. Ale w pesymistycznym
przypadku mnożenie zje ci pięć cyfr znaczących. Teraz wystarczy że od
obu wyników odejmiesz 1.0001 i zostaną ci tylko okruszki.
-
63. Data: 2014-08-05 09:11:29
Temat: Re: szybki logarytm
Od: Tomasz Kaczanowski <kaczus@dowyciecia_poczta.onet.pl>
W dniu 2014-08-05 08:38, slawek pisze:
> On Tue, 05 Aug 2014 08:00:26 +0200, Tomasz Kaczanowski
> <kaczus@dowyciecia_poczta.onet.pl> wrote:
>> Mnożenie jest operacją szybką i w zasadzie jedyna możliwość utraty
>> dokładności to konwersja inta na doubla (i na odwrót) - więc to
> zależy
>
> Oblicz 1.000000001 do potęgi 100000 przez kolejne mnożenia i wzorem.
> Oczywiście przez mnożeniem będzie szybciej jeżeli masz stabilicować
> wartości dla wykładników od 1 do 100000. Ale w pesymistycznym przypadku
> mnożenie zje ci pięć cyfr znaczących. Teraz wystarczy że od obu wyników
> odejmiesz 1.0001 i zostaną ci tylko okruszki.
wypowiadałem się tylko co do samych mnożeń - odejmowanie i dodawanie
wprowadzają duże zniekształcenia - ze względu na sprowadzanie danych do
wspólnego wykładnika. Druga rzecz, to zaokrąglenia, jeśli skończy się
zakres, ale tego nie przeskoczy się również w intach :) Jak napisałem -
zależy do czego się potrzebuje wyników, taka metoda będzie najlepsza.
--
Kaczus
http://kaczus.ppa.pl
-
64. Data: 2014-08-05 21:51:57
Temat: Re: szybki logarytm
Od: "slawek" <h...@s...pl>
Użytkownik "Tomasz Kaczanowski" napisał w wiadomości grup
dyskusyjnych:53e083a0$0$2146$6...@n...neostrada
.pl...
>zakres, ale tego nie przeskoczy się również w intach :) Jak napisałem -
>zależy do czego się potrzebuje wyników, taka metoda będzie najlepsza.
IMO, trudno jest zrobić "lepszy logarytm" niż jest on w hardware FPU. Chyba
że nie ma FPU, albo mamy jakieś szczególne fochy.
Owszem, opłaca się np. stablicować logarytm... jeżeli np. wiemy że dane
wejściowe to liczby od 1.00 do 2.00 z krokiem 0.01. Bo "ogólna funkcja
logarytm" nijak nie wie o tym, że nie ma być ogólna. Nie potrafi zrozumieć
czego my chcemy. Wiec wtedy opłaca się coś samemu robić - może akurat będzie
lepiej dopasowane do wymagań.
Gdzieś na półce w bibliotece (ale w wakacje jest zamknięta) leży świetna
książeczka z lat 50-tych, której autorzy zachwycają się pracą jakiegoś
Anglika (chyba lata 30, może 40), który stablicował współczynniki
wielomianów niskiego stopnia dających "całkiem dobre przybliżenia" funkcji
elementarnych. Myślę, że po użyciu log(x) = log(y) + log(10^n) aby przejść
od dowolnego x do takiego y które należy do przedziału (1,10), a potem
podobna sztuczka aby mieć (5,10), można bez trudu np. znaleźć wielomian
Czebyszewa dostatecznie dobry aby udawał logarytm.
AFAIR, w biblioteczce CEPHES była i funkcja logarytm:
http://www.boutell.com/lsm/lsmbyid.cgi/000626
Intel oferował jakąś "szybką bibliotekę" Approximate Math Library:
https://software.intel.com/en-us/articles/avoid-bott
lenecks-in-simple-math-functions
Z drugiej strony, wbudowany w FPU logarytm i tak wystarcza w 99% (lub
bardziej), a szkoda tracić czas... który lepiej przyda się do czegoś innego
niż pisanie własnej implementacji logarytmu.
Problem z podziałką gdzieś kiedyś już mi wyskoczył: mocno się zdziwiłem
(swego czasu) że iterowanie x := x + dx potrafi się przecudnie rozjechać (na
double'ach!) dla nietypowych danych. Np. gdy się chce pokazać ma wykresie
jak wynik Czegoś-Tam asymptotycznie zbliża się do dokładnego: masz oś od
xmin do xmax, ale różnica pomiędzy nimi to np. 3*epsilon maszynowy.
Nota bene, miałem zawsze trudności z tabliczką mnożenia (mnożenie w słupkach
ma koszt O(N*M)), więc nauczyłem się na pamięć mantys logartymów: zamiast
mnożyć 2*5 dodajesz 3010 do 6990 i wychodzi 1. Po pewnym czasie nie chce ci
się już patrzeć na tablice (Wojtowicza), bo i tak wiesz co tam jest.
Kalkulator? No cóż, były arytmometry mechaniczne. Fajne, na korbkę, taki
steam-punk jakby dziś na to popatrzeć.
-
65. Data: 2014-08-09 20:13:24
Temat: Re: szybki logarytm
Od: bartekltg <b...@g...com>
W dniu wtorek, 5 sierpnia 2014 09:11:29 UTC+2 użytkownik Tomasz Kaczanowski napisał:
> wypowiadałem się tylko co do samych mnożeń - odejmowanie i dodawanie
> wprowadzają duże zniekształcenia - ze względu na sprowadzanie danych do
> wspólnego wykładnika. Druga rzecz, to zaokrąglenia, jeśli skończy się
> zakres, ale tego nie przeskoczy się również w intach :) Jak napisałem -
> zależy do czego się potrzebuje wyników, taka metoda będzie najlepsza.
>
Ale on tylko tradycyjnie trolluje.
Mnożenie 100 000 razy daje na doublach błąd 27x epsylon, użycie pow
(w tej implementacji, czyli standardowy 32bitowy mingwin) 2x eps.
Straszne to;-) Różnica taka potworna, a w końcu generując punktu
do wykresu potrzebujemy ich 100 000...
A przynajmniej straszne się wydaje, póki nie powiem, że wynik
porównywałem nie z wynikiem uzyskanym z oryginalnej liczby,
a z liczby jaka w double reprezentuje "1.000000001".
Jej wartość to 1.000000001000000082740370999. Różnią się o < epsylon.
Ale jeśli porównam wyniki programu nie z
1.000000001000000082740370999^100000 a z dokładnym
1.000000001^100000 błąd wyniesie ponad... 37000 epsylonów.
Straszne i dziwne...
Nie takie straszne, bo to nadal tylko 10^-12.
Nie takie tez dziwne, bo pochodną x^n jest n*x^(n-1)
U nas w okolicach x=1 zmiana o dx zwiększa wartość o 100000dx,
błąd dyskretyzacji został rozdmuchany przez naszą funkcję
o kilka rzędów wielkości, mimo bardzo dobrego liczenia już
w arytmetyce fl.
pzdr
bartekltg
[znów z googla, sorry:)]
-
66. Data: 2014-08-09 20:18:50
Temat: [OT] szybki logarytm
Od: feldmarszałek tusk <N...@g...pl>
OT, jak zrobiłeś ten indeks górny?
-
67. Data: 2014-08-09 21:53:32
Temat: Re: szybki logarytm
Od: slawek <f...@f...com>
On Sat, 9 Aug 2014 11:13:24 -0700 (PDT), bartekltg
<b...@g...com> wrote:
> Mnożenie 100 000 razy daje na doublach błąd 27x epsylon, użycie pow=
Jak to oszacowaleś?
-
68. Data: 2014-08-09 22:09:39
Temat: Re: szybki logarytm
Od: slawek <f...@f...com>
On Sat, 9 Aug 2014 11:13:24 -0700 (PDT), bartekltg
<b...@g...com> wrote:
> do wykresu potrzebujemy ich 100 000...
Zależy kto co robi, np. równania rózniczkowe cząstkowe: kilka
miliardów operacji to raczej niedużo, to tylko parę sekund pracy CPU.
-
69. Data: 2014-08-11 10:39:50
Temat: Re: szybki logarytm
Od: "slawek" <h...@s...pl>
Użytkownik "slawek" napisał w wiadomości grup
dyskusyjnych:a...@n...v.pl
...
>On Sat, 9 Aug 2014 11:13:24 -0700 (PDT), bartekltg <b...@g...com>
>wrote:
>> Mnożenie 100 000 razy daje na doublach błąd 27x epsylon, użycie pow=
>Jak to oszacowaleś?
Tzn. nie że tyle ci akurat wyszło dla jakiś danych (czemu nie? mogło).
Ale że:
sto_tysiecy = 100000;
double a[1..sto_tysiecy] := jakies_random(&seed);
result := jakies_random(&seed);
for ( i in 2..sto_tysiecy ) loop
result := result * a[i];
end
error := 27*epsilon; // see: Bartekltg Great Theorem
Jak uzasadnisz swoją regułę, tj. ostatnią linijkę w pseudokodzie wyżej
zwłaszcza gdy np. "przypadkiem" liczby a[n] := n/1000 ???
IMHO, jeżeli a[n] jest bliskie 1, np. dla każdego n = 1,2, ..., m mamy 0.99
< a[n] < 1.01, to w pesymistycznym przypadku można spodziewać się błędu
m*epsilon.
Dlaczego tak? Bo przy mnożeniu "z grubsza" dodają się błędy względne, a te
dla a[n] prawie równej 1 są prawie równe błędom bezwzględnym. Dodawanie jest
w sensie abs-ów, bo - w pesymistycznym przypadku - wszystko będzie przeciwko
tobie. (W optymistycznym przypadku błędy się całkowicie wyzerują. Tak na
czuja przeciętnie powinno wyjść jakieś sqrt(m)*epsilon, co przy m równym 100
tysięcy daje 300*epsilon.)
Jak jeszcze tego nie rozumiesz, to prosty przykład: a = 17.0, b = 19.0, c =
23.0E+300. Załóżmy że wynik mnożenia a*b, obarczony błędem epsilon,
zapisujesz w rejestrze R. Załóżmy, że mnożenie R*c jest dokładne (jeżeli nie
jest, to jest jeszcze gorzej, ale optymistyczny przypadek). Wtedy błąd
mnożenia a*b powiększa się c krotnie, czyli wynosi epsilon*23.0E+300.
Ponieważ epsilon jest rzędu 1.0E-16, to masz błąd końcowego wyniku ponad
1.0E280. Ździebko dużo w porównaniu do 27*epsilon.
Względny błąd jest ok, ale reguła 27*epsilon po prostu nie działa.
Tak samo jak nie działa reguła "liczby nieparzyste są pierwsze" - i nie
pomaga że 13 jest nieparzysta i pierwsza.
-
70. Data: 2014-08-11 10:49:29
Temat: Re: szybki logarytm
Od: "slawek" <h...@s...pl>
Użytkownik "Borneq" napisał w wiadomości grup
dyskusyjnych:lr7uc5$oae$...@n...news.atman.pl...
>Skąd wziąłeś te wzory? Czy podobnie daje się wyprowadzić na inne funkcje? A
>jeśli chodzi o wolniejsze szeregi ale pozwalające wyliczyć z dowolną
>dokładnością, to gdzie można natknąć się na zbiór takich wzorów? Trochę
>jest na exp, sinus,czy cosinus w tablicach matematycznych.
Obejrzyj sobie:
http://www.springer.com/mathematics/applications/boo
k/978-0-387-48806-6
Za friko jest (w domenie publicznej):
http://people.math.sfu.ca/~cbm/aands/
Czasem przydaje się:
http://functions.wolfram.com