better floating point printing routines

kris gave me a nice paper about conversion to base 10 it with code in scheme, I will try to convert it...

I have converted the core of the printing algorithm of

www.cs.indiana.edu/~burger/FP-Printing-PLDI96.pdf

in D (not optimized, with focus on correctness).

To use this I would need the following things from BigInt? * number of digits (there is numBytes, but said to be deprecated and I suspect does not remove the leading zeros. * division and reminder (as this are calculated together optimally one would have a combined in an unique function, not in two calls). * simple conversion to from ints

If we focus just on base 10 (B=10) we can get rid of the duplicated factors 2, thereby reducing the size of the numbers.

So from the geekiness point of view it would be a very nice addition, I am not aware of anybody (aside from scheme) that really does this, seeing 0.78 printed as 0.78 and not or something like 0.78000000000000003 would be great. On the other hand this algorithm uses arbitrary integers (that Don implemented), but these integer might get rather large for reals, as the exponent might be quite large/small, and so the number involved. So I fear a bit the outcome

What do you think?

Just to give an idea of the number involved in the worst case: the floating point number is represented as fraction, with the smallest/biggest exponent the numerator/denominator can get quite big, ignoring the mantissa (which adds 1-2 or at most 4 bytes), and +/-1

float: emax=2**7-1=127 maxbytes no reduction ~ (2**7)/4=32 maxbytes reduction ~ (2**7)*(1-log(2)/log(10))/4=23

double emax=2**10-1=1023 maxbytes no reduction ~ (2**10)/4=256 maxbytes reduction ~ (2**10)*(1-log(2)/log(10))/4=179

80 bit emax=2**14-1=16383 maxbytes no reduction ~ (2**14)/4=4096 maxbytes reduction ~ (2**14)*(1-log(2)/log(10))/4=2863

quad precision emax=2**14-1=16383 maxbytes no reduction ~ (2**14)/4=4096 maxbytes reduction ~ (2**14)*(1-log(2)/log(10))/4=2863

So this are the numbers, Don can you maybe comment how feasible it would be to handle (in the extreme cases of 80bit and quad precision integer numbers (assuming reduction) of 2863, i.e. 2863*4+mantissa=11452+(64 or 113) bits? For doubles, I think this is feasible.

What do you mean by 'feasible'? Are you worried about memory allocation, or speed? Which operations are involved? Also, when you say you need 'number of digits' with leading zeros, what do you mean? I'd like to change the interface to allow you to implement this.

BTW, one of the reasons I started on BigInt? was for use in floating point output. Then I got a bit more interested in it in general.

with feasible I mean both computational cost and allocation, allocation I think is rather large for just converting, but is just a few KB.

I found this article on the other direction (reading numbers), that says that the routines by Gay (that gave up arbitrary arithmetic, I think) are normally much faster, and almost always correct.

ftp://ftp.ccs.neu.edu/pub/people/will/retrospective.pdf

I haven't looked in detail and we are interested in the other direction, but maybe there are obvious improvements to this algorithm (apart getting rid of the common powers of 2).

With number of number of digits I mean the number of binary digits to the leading 1 (not the number of bytes that would include the leading zeros.

ok I checked a little better the things. Basically Gay routines have better heuristic, and use floating point arithmetic when possible. This mekes them very fat if just few digits are requested. If one wants all the digits then the algorithm I translated is better, and is always exact.

I would like to have something like that, I think that the better heuristic of Gay for k can be used, so we would have an algorithm that would work for basically all floating point types, and would print the exact number of relevant digits in any case. Having 0.78 printed as (I am just inventing the number of zeros) 0.7800, 0.7800000, 0.780000000 depending if it is a float, double or real is (for me) deeply satisfying, and if we do it, it should become the default printing format.

Again what we would need is: - int (or long) to BigInt?, the nicest way would be to have BigInt?(i), and BI<<32+BigInt?(i), and the other direction i=BI&0xFFFF_FFFF; BI>>=32; - fast multiplication *10 (if it is worth it) - division of two large BigInts?, giving a reminder (large), and the quotient (a 0-9 digit).

Fawzi: early version of Tango did use the library from David Gay for FP input and output, though it was abandoned for a number of reasons. It doesn't support 80bit, which was problematic for some.

I think that at least some of the heuristic should be used for a small number of digits, and the general routine would use the algorithm by Burger and Dybvig. The only part of the algorithm that I am not sure if valid for 80bit and quad is the guess for k, but Gay has a better guess for it, so using that I think we are safe also for them (it should be checked, but I believe it will be so).

It seems that TCL did/wanted to do something similar (but only for doubles)

http://www.tcl.tk/cgi-bin/tct/tip/132.html

yeah, that makes sense. Interesting TCL link

core of the correct printing algorithm (with english labels, and scheme primer)

The low-level functions (eg in math.internal.BiguintCore?) will probably be more useful in the end than BigInt? itself, since they don't use memory allocation and have simple algorithms (and they give you div and mod simultaneously).

In terms of performance, multiplication n bytes*m bytes is about 6*n*m clock cycles for small m. Division n bytes/m bytes is in ~50*m*(n-m). It's not so bad if m is almost as large as n, or if it's small.

A bigger concern to me is how big the table of powers needs to be. For floats, you could have an entry for each of the 256 possible exponents, but for 80-bit reals the table gets huge.

There are a couple of factors which might make this whole thing more feasible: (1) D uses IEEE arithmetic. Consequently we have the nextUp/nextDown functions which TCL, scheme, etc don't have. (2) It _might_ make sense to have the lowest-level BigInt? functions (add, subtract, multiply-and-add) in core -- provided you can prevent them getting linked into every 'hello world' app. They actually aren't terribly long once compiled, and it IS something which you'd want to redo for every architecture. There's a fair bit of very useful stuff you can do with them.

David Gay's code is partly based on the fact that fast BigInt? functions are not available. I'm pretty intrigued by this whole thing, since I've never heard of anyone doing correct 80-bit printing before.

Yes it makes definitely sense to use BiguintCore?.

I was discussing with kris, I think that it makes sense to make the FP formatter pluggable. The code can be written using the type U, that can be a real of BigInt?, and then one can have converters using floats, BigInts?, and a hybrid switching between both (using Gay's strategy).

I was wondering if it was not possible to find a fast way to calculate the bit pattern of 5^{k (which is what we really need), but in the worst cast it is just log(k) to calculate that... so I would tabulate just the few closest...}

I have seen nextUp/nextDown, and thought about using them to find v+/v-.

... but I do not know if it is really useful. (was missing from the previous message

any further progress on this one?