Karim.Belabas on Wed, 4 Apr 2001 11:16:41 +0200 (MET DST) |
[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]
Re: [GP/PARI-2.1.0] arithmetic weirdness |
On Wed, 4 Apr 2001, Gerhard Niklasch wrote: [...] > The exact result of the multiplication %47*100 would have five extra > 1 bits to the right, which we're pushing out during normalization. > > We shouldn't drop them, we should remember enough of them to round > upward. The standard guard bits + sticky bit technique is being > called for here. > > For elementary operations and for as many algebraic and transcen- > dental functions as feasible, we should return the correctly rounded > t_REAL representation, at the appropriate precision, of the exact > result of applying the operation (in the reals, not in the t_REALs) > to the input value(s). Perhaps stuff like LLL would be a trifle > better behaved if we didn't drop bits in simple multiplications... > It is not necessary to compute the full 192-bit product of two 96-bit > significands to achieve this; a few extra bits suffice to do the > tracking. This is more or less what is being done internally. All transcendental functions are computed with one extra word of precision, then the result is truncated. The first two things (t_REAL) LLL does are 1) compute the number of significant words used in the input 2) add one guard word to all entries [ x = gmul(x, realun(prec+1)) ] It's a(n arguable) design choice, for the sake of maximal speed for low precision inputs: the t_REAL type corresponds to low-level multiplication; correct rounding and precision settings are the caller's responsibility. Of course you may end up with a number of recursive calls adding up their own extra word of precision until sluggish death occurs [ remember the zetak soaring accuracy problem ? ], so this places another burden on the PARI programmer [ stack management being the number 1 annoyance, until you get really used to it. Ubiquitous typecasts as opposed to, eg. structs or OO approach being my number 2 ]. It would be nicer if PARI also supported directed rounding or interval arithmetic; very useful in a number of settings. Such a package was submitted 3 years ago by Manfred Radimersky but I don't like changing the kernel [e.g Karatsuba multiplication for t_REAL has been there for 4 years, and is still not enabled by default, nor is it documented!] and, to my great regret, I unfortunately never had to time to check and adapt/include it. It still is on my personal TODO list. Karim. -- Karim Belabas email: Karim.Belabas@math.u-psud.fr Dep. de Mathematiques, Bat. 425 Universite Paris-Sud Tel: (00 33) 1 69 15 57 48 F-91405 Orsay (France) Fax: (00 33) 1 69 15 60 19 -- PARI/GP Home Page: http://www.parigp-home.de/