Karim BELABAS on Tue, 25 Feb 2003 00:21:28 +0100 (MET)

[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]

 Re: zetakinit() puzzle

```On Mon, 24 Feb 2003, Igor Schein wrote:
> \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\
> ? plot(x=2,3,zetak(zetakinit(y),x))
>
> 1.6449341 """""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""`
[...]
> -2.39e+18 |.............................................................._
>           2                                                              3
> \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\
>
> I expect the graphs to be identical.

Catastrophic precision loss [ due to the computation of Gamma(1 - s), for s
very close to an integer ]:

(23:56) gp > zetak(zetakinit(y), 3-1e-28)
%1 = -719751868412563322355.4757780

I have fixed an obvious mistake in the computation (and the above symptom),
but there is still a huge instability if zetakinit() was not computed to a
sufficient accuracy:

(23:56) gp > \p28
(23:56) gp > zetak(zetakinit(y), 3-1e-50)
%2 = -8.97164492 E34 \\ bogus
(23:56) gp > \p50
realprecision = 57 significant digits (50 digits displayed)
(23:57) gp > zetak(zetakinit(y), 3-1e-50)
%3 = 1.2020569031595942854001807100992369873297281694155
\\ better, but only the first 20 or so digits are correct

This function should be rewritten. But, as far as I understand the method,
precomputing zetakinit independently of the points where zeta_K is to be
computed simply doesn't work. It will give sensible results when applied to
sensible regions [ where 'sensible region' is defined by the fact you get
consistent results when increasing the precision ]. There are proven error
bounds, but they depend on s, and PARI can't apply them in zetakinit()...

Cheers,

Karim.
--
Karim Belabas                    Tel: (+33) (0)1 69 15 57 48
Dép. de Mathématiques, Bât. 425  Fax: (+33) (0)1 69 15 60 19
Université Paris-Sud             Email: Karim.Belabas@math.u-psud.fr
F-91405 Orsay (France)           http://www.math.u-psud.fr/~belabas/
--