Re: real polynomials

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

* John Cremona [2015-11-06 16:38]:
[...]
> I will give you a laugh by showing you my own version:
> 
> GEN dbltorat(double x) // assumes x=a/2^e
> {
>   long e;
>   GEN res = cgetg(3, t_FRAC);
>   if (x==0)
>     {
>       gel(res,1) = cgeti(0);
>       gel(res,2) = cgeti(1);
>     }
>   else
>     {
>       gel(res,1) = mantissa_real(dbltor(x), &e);
>       gel(res,2) = cgeti(2<<e);
>     }
>   return res;
> }
> 
> (tested!  a little)

There are a few problems indeed:

1) cgeti(0) / cgeti(1) / cgeti(2<<e) don't do what you think

In fact, both 0 and 1 are illegal arguments corrupting the stack
(immediately for 0, after other low-level commands such as setsigne() for 1) !


stoi(0), stoi(1) will both work [ think signed_long-to-integer ].
And so would utoi(0), utoi(1)   [ think unsigned_long-to-integer ]

cgeti(n) is a low-level routine reserving n words for later assignment
of an integer, it does not return a valid GEN in itself (and we must
have n >= 3 ) !

2) you probably meant  1 << e instead of 2 << e;  

3) int2n(e) will be more robust than stoi(1 << e), which breaks for 
  e >= BITS_IN_LONG-1 

4) t_FRAC members must be in normal form (lowest terms). Just dividing
the numerator by the denominator will do that automatically:

  num = mantissa_real(dbltor(x), &e);
  den = int2n(e);
  return gdiv(num, den);



Use high level constructors (such as gdiv / stoi) as much as possible, instead
of low lever allocations (such as cgetg / cgeti); the latter is marginally more
efficient but leads to (many) pitfalls...

Don't bother with garbage collecting for first programs; if absolutely
necessary (stack overflows), stick to the

  pari_sp av = avma;
  ...
  return gerepilecopy(av, result);

idiom. Good luck !

Cheers,

    K.B.

P.S. There's a PARI/GP Atelier in Grenoble from January 11th-15th 2016.
It's the ideal event to quickly learn the basics :-)

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