W dniu 2013-03-26 12:25, M.M. pisze:
> Niczego to nie zmienia. Logarytm tez mozna policzyc na mantysie o 2 bity
> dluzszej (np. na precyzji 66bitow) i wynik bylby dokladny.


:)

Szkolny przykład.
x=1/4;
x <- 4*x*(1-x).

Powinien wyjść cykl 1/4, 3/4, 1/4...


for (int i=0;i<1000;i++)
x = 4*x*(1-x);

#include <cstdio>
#include <cmath>
#include <time.h>
#include <limits>

double bla(double x,const int n)
{
for (int i=0;i<n;i++)
{
x= 4*x*(1-x);
}
return x;
}

int main()
{
double a1 = 1.0/4.0;
double a2 = a1+a1*std::numeric_limits<double>::epsilon();
double b1,b2;

for (int n=0; n<=70;n+=10)
{
b1= bla(a1,n); b2= bla(a2,n);
printf("%d %lf %lf %lg\n", n, b1, b2, (b1-b2) );
}
return 0;
}


0 0.250000 0.250000 -5.55112e-017
10 0.750000 0.750000 1.13687e-013
20 0.750000 0.750000 1.16415e-010
30 0.750000 0.750000 1.19209e-007
40 0.750000 0.749878 0.00012208
50 0.750000 0.616384 0.133616
60 0.750000 0.613975 0.136025
70 0.750000 0.995892 -0.245892

Ani dwa, ani nawet 2000 bitów nic nie pomoże.
Najmniejsza odchyłka i trajektorie geometrycznie
się rozbiegają.

Dla lepszego związku z tematem dodam, że dwa bity więcej
nie spowodują, że zawsze dostaniesz dobry wynik. Dwa
bity więcej oznaczają, że najczęściej dostaniesz poprawny
wynik.

Któryś logarytm którejś liczby dobrze zapisanej
w systemie dwójkowym zapisze się źle mimo tych
dodatkowych bitów i katastrofa gotowa.


Zmieńmy naszą iterację na
x= exp(log(4.0)+log(x)+log(1-x));

0 0.250000 0.250000 -5.55112e-017
10 0.750000 0.750000 8.52651e-014
20 0.750000 0.750000 8.73115e-011
30 0.750000 0.750000 8.9407e-008
40 0.749969 0.749878 9.15621e-005
50 0.718126 0.616384 0.101742
60 0.950183 0.613966 0.336216
70 0.840387 0.994657 -0.15427

Odpłynęliśmy równie daleko, mimo, że zarówno argument
naszego wyrażenia, jak i wynik były zapisywalne dokładnie.

pzdr
bartekltg