Re: Is There a Way to Rationalize a Decimal in Pari/GP?

[Kurt Foster <drsardonicus@earthlink.net>]

On Jun 5, 2009, at 8:07 AM, Rick Regan wrote:

> I didn't mention that I am only looking at this from a limited point  
> of view. I want to convert dyadic decimal fractions -- those that  
> represent double-precision floating-point binary values.

This suggests the following problem: Assuming the input x can be  
represented as a normalized double-precision value, find that value  
(or one such value).  Then, determine the simplest fraction which is  
consistent with that value.

Since a normalized double-precision value has 52 significant bits with  
an implied leading 1, the default precision (either 28 or 38 decimal  
digits) should be plenty to determine a normalized double-precision  
representation of x, assuming it has (at least) one.

Assuming x!=0, e can let k = 52 - floor(log(x)/log(2)), then N =  
round(x*2^k) satisfies 2^52 < x <= 2^53.  Assuming N < 2^53, the  
problem then becomes that of determining the simplest fraction which  
lies in the open interval

](N - 1/2)*2^(-k), (N + 1/2)*2^(-k)[

There are many details to deal with, of course.  Among other things,  
if x is such that x*2^k is exactly midway between two integers, N  
isn't uniquely determined.  Or what of N = 2^53?  Or x == 0? And, of  
course, there's checking that the exponent k is consistent with a  
normalized double-precision value.  But apart from these annoyances,  
there will be a unique answer.  There's a way to use simple continued  
fractions to find it, but (as with the double-precision  
representation) I'm too lazy to work out the details.

I did use  the above formulation on some examples.  The sample input  
from the thread

0.1000000000000000055511151231257827021181583404541015626

gives the fraction 1/10,

0.33333333333333333333333333333333333333333

gives the fraction 1/3, and

0.759765625

gives the fraction 389/512.