Karim Belabas on Thu, 02 Mar 2006 22:04:22 +0100


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

Re: qflllgram


* Gonzalo Tornaria [2006-02-27 18:47]:
> Strange behaviour:
> 
> ? qflllgram([1,0;0,-1],0)
>   *** qflllgram: the PARI stack overflows !
>   current stack size: 128000000 (122.070 Mbytes)
>   [hint] you can increase GP stack with allocatemem()
> 
> ? allocatemem()
>   *** allocatemem: Warning: doubling stack size; new stack = 256000000
> (244.141 Mbytes).
> ? qflllgram([1,0;0,-1],0)
>   *** qflllgram: length (lg) overflow

Junk in, junk out: the function is undefined on indefinite matrices.

More precisely:

For flag = 0, the implementation uses floating point computations and
increases the internal precision when incremental Gram-Schmidt fails
(basis vector with non-positive square norm).

When the input is exact, we enter an infinite loop, and eventually
overflow the capacities of the implementation. When the input is
inexact, we give up relatively quickly.

> It /does/ work with flag=1 (lllgramint):
> 
> ? qflllgram([1,0;0,-1],1)
>   *** qflllgram: not a definite matrix in lllgram

For flag = 1, we use the atrociously slow all-integer algorithm. But a
negative vector length in Gram-Schmidt is immediately flagged as an error:
since all computations are exact, it can't be due to roundoff errors.


I don't see why we should clutter the code and slow down the floating
point code to check the input in case the user would be playing tricks
on us. Try it at \g4 and a non trivial computation, like

  \g4
  \p300
  algdep(3^(1/5)+2^(1/6), 30)

Quite a number of precision changes, heh ? [ This could be performed
much more efficiently, but I haven't had time to implement Phong &
Stehle's L^2 algorithm yet ]

Cheers,

    K.B.
-- 
Karim Belabas                  Tel: (+33) (0)5 40 00 26 17
Universite Bordeaux 1          Fax: (+33) (0)5 40 00 69 50
351, cours de la Liberation    http://www.math.u-bordeaux.fr/~belabas/
F-33405 Talence (France)       http://pari.math.u-bordeaux.fr/  [PARI/GP]