Дисксуссия получила развитие
Цитата:
actually, the only command that increases the precision there is
the addition s + i; we add to the t_REAL s the t_INT i which has
- infinite accuracy (obviously)
- larger exponent than s after the first few loops
So the result s + i has a larger precision than s (in fact, the
bitprecision increases by 64 which is the minimal amount), raising that
to the power 1/5 keeps that precision. But the next loop with the new s
increases it again.
In practice, the precision increases by 64 bits every loop.
Cheers,
K.B.
K.B. это Karim Belabas, один из разработчиков.
Он говорит следующее. Поскольку целое число, представленное с бесконечной точностью
i складывается с плавающим
s, и при том
, то точность результата на всякий случай увеличивается на 64 бита (почему на 64 - потому что это кратность увеличения точности, т.е. минимальный шаг увеличения точности). Это происходит каждый цикл, соответствено, если циклов как в моём примере, который я им отправил, составляет 2000, то точность представления
s увеличивается в сумме на 64*2000=128000 бит. Ну это примерно 38 тыс десятичных цифр.
Ну что ж, проверяем:
Код:
? s=0.0;print(precision(s));n=20;for(i=1,n,s=(s+i)^(1/5);print("i=",i," precision=",precision(s)))
38
i=1 precision=38
i=2 precision=38
i=3 precision=38
i=4 precision=38
i=5 precision=38
i=6 precision=38
i=7 precision=38
i=8 precision=38
i=9 precision=38
i=10 precision=38
i=11 precision=38
i=12 precision=38
i=13 precision=38
i=14 precision=38
i=15 precision=38
i=16 precision=57
i=17 precision=77
i=18 precision=96
i=19 precision=115
i=20 precision=134
cpu time = 4 ms, real time = 4 ms.
?
Видим что действительно, начиная с
i=16, точность
s с каждым циклом увеличивается 19-20 дсятичных значащих цифр.
Это конечно удивительно на первый взгляд. Но далее Bill Allombert (тоже разработчик) приводит такой пример:
Цитата:
Alas, PARI t_REAL precision is always a multiple of 64 bits.
Whether this is pessimistic or optimistic is a matter of perspective:
Compare this:
? s=0.;n=20;t=100.;for(i=1,n,s=(s+t)^(1/5));forstep(i=n,1,-1,s=s^5-t);s
%38 = -3.3995948126355403614265653237903756818E64
? s=0.;n=20;t=100;for(i=1,n,s=(s+t)^(1/5));forstep(i=n,1,-1,s=s^5-t);s
%39 = -4.3414201952566371428595921256072201009E-56
%38 is completly wrong while %39 is quite accurate.
Cheers,
Bill
Что мы тут видим: сперва считаем сумму, потом пытаемся всё вычесть, должны получить ноль.
Если не увеличивать точность (инкремент в основании делаем плавающее 100.0), то в результате получаем
что от нуля, мягко говоря, далековато. А если точность увеличивать (инкремент в основании делаем целочисленный 100), то получаем
что вполне к нулю близко. Соотношение результатов -- более 100 десятичных порядков.
По крайней мере, причина замедления и поедания памяти теперь понятна, ну а каждый решает стоило ли оно того в его конретном случае.
Какой из этого вывод? Ну... вывод простой - мойте руки перед едой, соблюдайте гигиену вычислений с плавающей точкой