Floating point exactness, another alternative.

Steven Pemberton steven at cwi.nl
Mon Aug 13 22:36:38 AEST 1990


In article <713 at tetrauk.UUCP> rick at tetrauk.UUCP (Rick Jones) writes:
> Thanks to all those who responded to my query about exactness (or lack of it)
> in FP numbers.  It has confirmed my suspicion that there isn't a
> clean solution to this problem.
> 
> The responses can be grouped into 3 broad categories:
> a.	"The effective precision is a function of the algorithm and
> 	initial data, not just the representation"
> b.	"Analyse the bit-pattern of the IEEE format"
> c.	The "constructive reals" solution
> d.	The "continued fractions" solution

I'm afraid I missed the original query, but I feel I should mention
what ABC does, especially since it does solve the problem of
representing amounts of money exactly (which was the problem domain
for the query).

ABC has exact numbers, which it represents internally as (unbounded)
rational numbers, reduced to lowest terms. So 1.25 which is (125/100)
is represented internally as (5, 4). You can always round to some
number of digits if you want. So "1/3" is (1, 3), and "2 round (1/3)"
is (33, 100). Arithmetic on exact numbers in general yields exact
numbers where it can, but functions like square root round to machine
accuracy.

For instance, here is a function (rep) that represents a number exactly,
indicating repeating decimals between square brackets (so 1/3 gives
0.[3]).

First an example of it in use:

	FOR i IN {1..20}: WRITE i, ": ", rep(1/i) /
	1 : 1
	2 : 0.5
	3 : 0.[3]
	4 : 0.25
	5 : 0.2
	6 : 0.1[6]
	7 : 0.[142857]
	8 : 0.125
	9 : 0.[1]
	10 : 0.1
	11 : 0.[09]
	12 : 0.08[3]
	13 : 0.[076923]
	14 : 0.0[714285]
	15 : 0.0[6]
	16 : 0.0625
	17 : 0.[0588235294117647]
	18 : 0.0[5]
	19 : 0.[052631578947368421]
	20 : 0.05

Here is the function:

	HOW TO RETURN rep x: \ The exact representation of a number
	    IF x = 0: RETURN "0"
	    PUT "`floor x`.", 10*(x mod 1), {} IN t, rem, index
	    WHILE rem not.in keys index:
		PUT t^"`floor rem`", 10*(rem mod 1), #t IN t, rem, index[rem]
	    IF rem = 0: RETURN t|index[rem]
	    RETURN "`t|index[rem]`[`t at index[rem]+1`]"

An explanation of parts of this:

	The algorithm works in the way you would do it by hand: it
	works out the next digit of the fraction and the remainder
	that is left, and for each remainder the position in the
	string where that remainder occurred. If the remainder has
	already occurred, then we have a recurring fraction, and we
	can stop.

	't' is the string we are building up, 'rem' the remainder we
	have now, and 'index' the record of where each remainder occurs.

	Within a string, parts in backquotes `` are evaluated.
	So "2**10 = `2**10`" gives "2**10 = 1024".

	The operator ^ joins two strings. The operators | and @ trim a
	string: t|4 returns a string consisting of the first 4
	characters of t; t at 5 returns the tail of t starting at the
	fifth character.

	x mod 1 gives the fractional part of x.

	The variable 'index' is a 'table', a generalisation of arrays.
	This one saves the position that each remainder occurs at.
		PUT #t IN index[rem]
	saves the length of the string t (#t) in the table at key 'rem'.

	'keys index' gives a list of the keys of the table, so
		WHILE rem not.in keys index:
	stops when we have a remainder that has already been handled.

Steven Pemberton, CWI, Amsterdam; steven at cwi.nl
"Let us go then you and I/while the night is laid out against the sky/like a
					smear of mustard on an old pork pie"
Nice poem Tom. I have ideas for changes though, why not come over? - Ezra



More information about the Comp.lang.c mailing list