Re: zetakinit() puzzle

[Karim BELABAS <Karim.Belabas@math.u-psud.fr>]

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/
--
PARI/GP Home Page: http://www.parigp-home.de/