about integer and float operations

Hi. During some experiments with Python and Haskell I found some important differences about how some integer and float operations are implemented. The first difference is about a `mod` b, when a and b are Float types. Python use the fmod function, and it also implement divmod; Haskell seems to lack support for this operation. The second difference is about the division of two integers. Consider this Python code:
from __future__ import division def fac(n): ... return reduce(lambda a, b: a * (b + 1), xrange(n), 1) ... float(fac(777)) Traceback (most recent call last): File "<stdin>", line 1, in <module> OverflowError: long int too large to convert to float fac(777) / fac(777) 1.0
Here CPython does not convert the two integers to float before to divide them, but make use of a special algorithm. GHC, instead, returns NaN I don't know if the implementations of divMod and "true" integer division in CPython are "robust", but there is some important reason why these operations are not supported in Haskell? Thanks Manlio Perillo

Manlio Perillo wrote:
The first difference is about a `mod` b, when a and b are Float types. Python use the fmod function, and it also implement divmod; Haskell seems to lack support for this operation.
Yes, Haskell does not implement the full IEEE. There are differing opinions about that: some say we should add it all in, and some say we should take it all out. Floating point is ugly, but useful. It seems that enough of IEEE is implemented for the vast majority of applications, so far.
fac(777) / fac(777) 1.0 Here CPython does not convert the two integers to float before to divide
The second difference is about the division of two integers. them, but make use of a special algorithm. GHC, instead, returns NaN
No, actually here Haskell shines. Perhaps this GHCi session will illuminate the issue for you: Prelude> let fac n = product [2..n] Prelude> fac 777 `div` fac 777 1 Prelude> fac 777 / fac 777 NaN In Haskell, the various integral and floating point types are completely separate. The literal 777 is overloaded - its type is 777 :: Num a => a so it can be used for both integral and floating point types. But the division operators are separate: div :: Integral a => a -> a -> a (/) :: Fractional a => a -> a -> a So when you use 777 together with /, the 777 is interpreted as a Fractional (defaulting to Double). And when you use it with div, the 777 is interpreted as an Integral (defaulting to Integer). Hope this helps, Yitz

Manlio Perillo wrote:
The first difference is about a `mod` b, when a and b are Float types. Python use the fmod function, and it also implement divmod; Haskell seems to lack support for this operation.
I wrote:
Yes, Haskell does not implement the full IEEE.
I spoke too soon. Data.Fixed.mod' supports this operation. I doubt that it is implemented via the IEEE primitive though. -Yitz

Yitzchak Gale ha scritto:
Manlio Perillo wrote: [...]
fac(777) / fac(777) 1.0 Here CPython does not convert the two integers to float before to divide
The second difference is about the division of two integers. them, but make use of a special algorithm. GHC, instead, returns NaN
No, actually here Haskell shines. Perhaps this GHCi session will illuminate the issue for you:
Prelude> let fac n = product [2..n] Prelude> fac 777 `div` fac 777 1 Prelude> fac 777 / fac 777 NaN
No, this is not as complete as it is done in Python. In recent versions of Python you have two division operators. The `/` operator *always* perform a true division. As an example, the division of two integers return a float. The `//` operator *always* perform a "floor" division. This happens for both integers and floats:
2.5 // 1.5 1.0 2.5 / 1.5 1.6666666666666667
In Haskell: Prelude> 2.5 `div` 1.5 <interactive>:1:0: Ambiguous type variable `t' in the constraints: `Integral t' arising from a use of `div' at <interactive>:1:0-12 `Fractional t' arising from the literal `1.5' at <interactive>:1:10-12 Probable fix: add a type signature that fixes these type variable(s) (but this seems to be available with Data.Fixed, however on Debian Etch I still have GHC 6.8.2). As for your example:
Prelude> let fac n = product [2..n] Prelude> fac 777 `div` (4 * fac 776) 194
This is incorrect, because `div` returns an integer, but I want a float with the exact result (194.25 with Python). If I'm correct, there is no operator/function, in Haskell, that perform an exact division between two integers and return a float: exactDiv :: (Integral a, Real b) => a -> a -> b I personally prefer the Python solution, where we have two operators with the same behaviour over all the numbers. In Haskell, something like (/) :: (Num a, Real b) => a -> a -> b (//) :: (Num a, Integral b) => a -> a -> b
[...]
Thanks Manlio Perillo

Manlio Perillo ha scritto:
[...]
I personally prefer the Python solution, where we have two operators with the same behaviour over all the numbers.
In Haskell, something like
(/) :: (Num a, Real b) => a -> a -> b
This should be (/) :: (Num a, Fractional b) => a -> a -> b but I'm not sure it is correct.
(//) :: (Num a, Integral b) => a -> a -> b
Manlio Perillo

Prelude> let i2fDiv a b = fromIntegral a / fromIntegral b Prelude> :t i2fDiv i2fDiv :: (Integral a, Fractional b, Integral a1) => a -> a1 -> b Prelude> 10 `i2fDiv` 3 3.3333333333333335 That what you're looking for? -Ross On Feb 4, 2009, at 4:22 PM, Manlio Perillo wrote:
Manlio Perillo ha scritto:
[...] I personally prefer the Python solution, where we have two operators with the same behaviour over all the numbers. In Haskell, something like (/) :: (Num a, Real b) => a -> a -> b
This should be (/) :: (Num a, Fractional b) => a -> a -> b
but I'm not sure it is correct.
(//) :: (Num a, Integral b) => a -> a -> b
Manlio Perillo _______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

Ross Mellgren ha scritto:
Prelude> let i2fDiv a b = fromIntegral a / fromIntegral b Prelude> :t i2fDiv i2fDiv :: (Integral a, Fractional b, Integral a1) => a -> a1 -> b Prelude> 10 `i2fDiv` 3 3.3333333333333335
Prelude> fac 777 `i2fDiv` fac 777 NaN
That what you're looking for?
No. I'm looking for an exact integer division that avoids overflows, if possible.
-Ross
Manlio Perillo

Manlio Perillo wrote:
No. I'm looking for...
Manlio - can you describe exactly what you want? Do you know exactly what you want? You have said that you want division like in Python - but that even that is not well defined: Python 2.6.1
3/5 0
Python 3.1
3/5 0.59999999999999998
Please tell us *exactly* what you want your division to do, on integers and on floating point. Or tell us you want division like / in Python 3.1, or whatever. -Yitz

Yitzchak Gale ha scritto:
Manlio Perillo wrote:
No. I'm looking for...
Manlio - can you describe exactly what you want? Do you know exactly what you want?
You have said that you want division like in Python - but that even that is not well defined:
Python 2.6.1
3/5 0
You have to: from __future__ import division In Python 3.x this is now the default behaviour
Python 3.1
3/5 0.59999999999999998
Please tell us *exactly* what you want your division to do, on integers and on floating point.
I'm posting more details in another response. Thanks Manlio Perillo

On 5 Feb 2009, at 10:38 am, Manlio Perillo wrote:
I'm looking for an exact integer division that avoids overflows, if possible.
What this sounds like to me is a request that the Prelude function 'fromRational' should work well. Since "The floating point literal f is equivalent to fromRational (n Ratio.% d), where fromRational is a method in class Fractional and Ratio.%constructs a rational from two integers, as defined in the Ratio library. The integers n and d are chosen so that n/d = f." If you cannot divide two Integers n, d accurately using (fromRational (n Ratio.% d) :: Double) that casts doubt on the trustworthiness of floating point literals. Suppose we have a function decodeIntegerAsFloat :: RealFloat a => Integer -> (Integer,a) such that if (s,m) = integer_to_scaled_float x then either x = 0 and s = 0 and m = 0 or x = m * 2**s (mathematically) and abs m \in [0.5,1.0). Then integer_ratio_as_float :: Floating a => Integer -> Integer -> a integer_ratio_as_float p q = (mp/mq)*(2.0^(sp-sq)) where (sp,mp) = decodeIntegerAsFloat p (sq,mq) = decodeIntegerAsFloat q You'd actually use scaleFloat; if the difference sp-sq is outside the range of Int the answer is going to be a signed zero or a signed infinity anyway. decodeIntegerAsFloat would sit very well in the RealFloat class alongside its model, decodeFloat. It has other uses. For example, you can use it to compute logarithms of Integers with much less worry about overflow.

Richard O'Keefe ha scritto:
On 5 Feb 2009, at 10:38 am, Manlio Perillo wrote:
I'm looking for an exact integer division that avoids overflows, if possible.
What this sounds like to me is a request that the Prelude function 'fromRational' should work well.
Just found that it actually works well, with Ratio Integer. There is a specialized function in GHC.Float that does the right job.
[...]
You'd actually use scaleFloat; if the difference sp-sq is outside the range of Int the answer is going to be a signed zero or a signed infinity anyway. decodeIntegerAsFloat would sit very well in the RealFloat class alongside its model, decodeFloat. It has other uses. For example, you can use it to compute logarithms of Integers with much less worry about overflow.
By the way, in GHC.Float there is a (private): -- Compute the (floor of the) log of i in base b. -- Simplest way would be just divide i by b until it's smaller then b, -- but that would be very slow! We are just slightly more clever. integerLogBase :: Integer -> Integer -> Int Regards Manlio Perillo

Manlio Perillo wrote:
By the way, in GHC.Float there is a (private): integerLogBase :: Integer -> Integer -> Int
Yes, I have needed this function many times. Too bad it is not exposed. In this case, though, we only need base 2. For that, it would be nice if we could just read it directly from the internal representation of the Integer. -Yitz

| Manlio Perillo wrote: | > By the way, in GHC.Float there is a (private): | > integerLogBase :: Integer -> Integer -> Int | | Yes, I have needed this function many times. | Too bad it is not exposed. So use the libraries process to propose exposing it! Simon

Manlio Perillo wrote:
I'm looking for an exact integer division that avoids overflows, if possible.
Richard O'Keefe wrote:
What this sounds like to me is a request that the Prelude function 'fromRational' should work well... If you cannot divide two Integers n, d accurately using (fromRational (n Ratio.% d) :: Double) that casts doubt on the trustworthiness of floating point literals.
No, that works fine: Prelude Data.Ratio> fromRational $ (3*10^1000)%(2*10^1000+1) 1.5
Suppose we have a function decodeIntegerAsFloat :: RealFloat a => Integer -> (Integer,a) such that if (s,m) = decodeIntegerAsFloat x then either x = 0 and s = 0 and m = 0 or x = m * 2**s (mathematically) and abs m \in [0.5,1.0).
Yes, that is what Manlio wants. Sometimes you need to divide two very large Integers with a floating point number as result, without the overhead of constructing a Rational from them.
Then integer_ratio_as_float :: Floating a => Integer -> Integer -> a integer_ratio_as_float p q = (mp/mq)*(2.0^(sp-sq)) where (sp,mp) = decodeIntegerAsFloat p (sq,mq) = decodeIntegerAsFloat q
You'd actually use scaleFloat; if the difference sp-sq is outside the range of Int the answer is going to be a signed zero or a signed infinity anyway.
You have just duplicated the CPython interpreter source code snippet that Manlio linked to.
decodeIntegerAsFloat would sit very well in the RealFloat class alongside its model, decodeFloat. It has other uses.
I agree, though I'm not sure it needs to be a method.
For example, you can use it to compute logarithms of Integers with much less worry about overflow.
Actually, efficient integer-valued logarithms for Integers are exactly what you need to implement decodeIntegerAsFloat. Best would be if we could just read that off from the internal representation of an Integer. Would you like to file a GHC issue for that? In the meantime, we could already expose the integer log function in a library using an efficient algorithm. Thanks, Yitz

Yitzchak Gale ha scritto:
[...]
Suppose we have a function decodeIntegerAsFloat :: RealFloat a => Integer -> (Integer,a) such that if (s,m) = decodeIntegerAsFloat x then either x = 0 and s = 0 and m = 0 or x = m * 2**s (mathematically) and abs m \in [0.5,1.0).
Yes, that is what Manlio wants. Sometimes you need to divide two very large Integers with a floating point number as result, without the overhead of constructing a Rational from them.
By the way: it is possible to use a private constructor (via some special GHC flag?). I would like to do a quick performance check using the existing fromRational specialization by constructing a Rational directly. I know that Haskell allows declaration hiding for program safety, but sometimes this can be a nuisance.
[...]
Thanks Manlio Perillo

| By the way: it is possible to use a private constructor (via some | special GHC flag?). | I would like to do a quick performance check using the existing | fromRational specialization by constructing a Rational directly. | | I know that Haskell allows declaration hiding for program safety, but | sometimes this can be a nuisance. I think you mean that if you have module M(f) where f = ... g = ... precompiled in a library, is it possible to call g? No, I'm afraid it isn't, at least in GHC. For example, 'g' may no longer exist... it may have been inlined into 'f'. Simon

On Wed, Feb 4, 2009 at 1:09 PM, Manlio Perillo
In Haskell, something like
(/) :: (Num a, Real b) => a -> a -> b
You probably want (Real a, Fractional b) => a -> a -> b. Int is an instance of Real... Real is the class of types that can be converted to Rational. Then we can define (/.) :: (Real a1, Real a2, Fractional a) => a1 -> a2 -> a x /. y = fromRational $ toRational x / toRational y whose type is more general than the one I gave above, but we can restrict it to that by changing the signature, if you like.
(//) :: (Num a, Integral b) => a -> a -> b
Again, Num is inappropriate, but we can define something similar: (//) :: (Integral b, Real a, Real a1) => a -> a1 -> b x // y = floor $ toRational x / toRational y Hope that helps, Max

Max Rabkin ha scritto:
[...]
Then we can define (/.) :: (Real a1, Real a2, Fractional a) => a1 -> a2 -> a x /. y = fromRational $ toRational x / toRational y [...]
(//) :: (Integral b, Real a, Real a1) => a -> a1 -> b x // y = floor $ toRational x / toRational y
Hope that helps, Max
Yes, thanks. However there is still a *big* problem: it is inefficient. Here is a Python version of the Chudnovsky algorithm [1] for computing Pi: http://paste.pocoo.org/show/102800/ On my system it takes 10 seconds. Here is an Haskell version: http://paste.pocoo.org/show/102801/ On my system it takes 30 seconds. Thanks Manlio Perillo [1] http://upload.wikimedia.org/math/6/6/8/6681cd21f3ca9bf13248a87d4202e06a.png

On Wed, Feb 4, 2009 at 2:03 PM, Manlio Perillo
Max Rabkin ha scritto:
[...]
Then we can define (/.) :: (Real a1, Real a2, Fractional a) => a1 -> a2 -> a x /. y = fromRational $ toRational x / toRational y
[...]
(//) :: (Integral b, Real a, Real a1) => a -> a1 -> b x // y = floor $ toRational x / toRational y
Hope that helps, Max
Yes, thanks.
However there is still a *big* problem: it is inefficient.
It sure is. To get the type you asked for, one *has* to go via rationals. But...
Here is an Haskell version: http://paste.pocoo.org/show/102801/
On my system it takes 30 seconds.
You're only dividing integers by integers to get Doubles. x /. y = fromIntegral x / fromIntegral y works just fine in this case. --Max

Max Rabkin ha scritto:
[...]
Here is an Haskell version: http://paste.pocoo.org/show/102801/
On my system it takes 30 seconds.
You're only dividing integers by integers to get Doubles.
x /. y = fromIntegral x / fromIntegral y
works just fine in this case.
No, this *does not works*. I get Nan as result, since fromIntegral a, when a is out of range of a double, returns Infinity. Prelude> (fromIntegral (fac 777)) :: Double Infinity CPython, does not try to convert to double, but, instead, to a scaled double.
--Max
Manlio Perillo

Manlio Perillo ha scritto:
[...] Here is an Haskell version: http://paste.pocoo.org/show/102801/
On my system it takes 30 seconds.
Sorry, I compiled without optimizations enabled. With -O2 it now runs in 20 seconds. What other optimizations should I enable? Thanks Manlio Perillo

Manlio Perillo wrote:
However there is still a *big* problem: it is inefficient.
Here is a Python version of the Chudnovsky algorithm [1] for computing Pi: http://paste.pocoo.org/show/102800/ On my system it takes 10 seconds. Here is an Haskell version: http://paste.pocoo.org/show/102801/ On my system it takes 30 seconds.
Ah, OK. Thanks. Now we have a well-defined problem. :) And that makes it clear that what you want is Python 3 division. Well, first of all, there are a few bugs in the Haskell version of your program - don't forget that unlike in Python, ranges in Haskell *include* the last number. Second, we should note that it is silly to compute 1000 terms in this sum. By the end, you are getting terms whose order of magnitude is 10 ^ (-15000). The Haskell version is spending a bit more time on those (useless) later terms. The reason is that in Haskell we are first reducing the fraction, then converting to Double. If you really did care about that amount of accuracy, you would leave the result as a Rational. Or as one of various other unlimited-precision floating point types that are available in Haskell libraries. This calculation in Haskell would take the same amount of time as it takes now. You would need to rewrite your Python program to do this, and I would guess it would run a lot slower. In our case, the Python division first does a quick estimate of the sizes of the two integers, and just returns zero if it sees that there will be underflow on conversion to double. So I made the following rough change to the Haskell: -- An exact division (/.) :: Integer -> Integer -> Double x /. y | y `div` x > 5*10^323 = 0 | otherwise = fromRational $ toRational x / toRational y Now the Haskell runs in about half the time as the Python on my machine. Obviously, the Haskell could easily be further optimized. -Yitz

Yitzchak Gale ha scritto:
Manlio Perillo wrote:
However there is still a *big* problem: it is inefficient.
Here is a Python version of the Chudnovsky algorithm [1] for computing Pi: http://paste.pocoo.org/show/102800/ On my system it takes 10 seconds. Here is an Haskell version: http://paste.pocoo.org/show/102801/ On my system it takes 30 seconds.
Ah, OK. Thanks. Now we have a well-defined problem. :)
Good :).
And that makes it clear that what you want is Python 3 division.
Well, first of all, there are a few bugs in the Haskell version of your program - don't forget that unlike in Python, ranges in Haskell *include* the last number.
Ah right, thanks!
Second, we should note that it is silly to compute 1000 terms in this sum. By the end, you are getting terms whose order of magnitude is 10 ^ (-15000).
Yes, of course. But I wrote the Python script as a response to a post on it.comp.lang.python, where an user was having overflow problems with a naive implementation of the algorithm. I have tested it with 1000 iterations just to check the raw performances. Then I tried to code it in Haskell, noting that there was no direct method for an exact division of two integer numbers.
The Haskell version is spending a bit more time on those (useless) later terms. The reason is that in Haskell we are first reducing the fraction, then converting to Double. If you really did care about that amount of accuracy, you would leave the result as a Rational.
This is not possible, since there is an irrational number: c ** (3. / 2.). Moreover the sum of rational numbers is rather expensive, in this case.
[...]
In our case, the Python division first does a quick estimate of the sizes of the two integers, and just returns zero if it sees that there will be underflow on conversion to double. So I made the following rough change to the Haskell:
-- An exact division (/.) :: Integer -> Integer -> Double x /. y | y `div` x > 5*10^323 = 0 | otherwise = fromRational $ toRational x / toRational y
Right, but I would like to see a proper implemented function for exact integer division in GHC.
Now the Haskell runs in about half the time as the Python on my machine. Obviously, the Haskell could easily be further optimized.
-Yitz
Thanks Manlio

On Thu, 2009-02-05 at 01:10 +0100, Manlio Perillo wrote:
Yitzchak Gale ha scritto:
In our case, the Python division first does a quick estimate of the sizes of the two integers, and just returns zero if it sees that there will be underflow on conversion to double. So I made the following rough change to the Haskell:
-- An exact division (/.) :: Integer -> Integer -> Double x /. y | y `div` x > 5*10^323 = 0 | otherwise = fromRational $ toRational x / toRational y
Right, but I would like to see a proper implemented function for exact integer division in GHC.
(%) *is* a proper function for exact integer division. But you'll find plenty of Haskellers to balk at calling anything that returns a Double `proper'. jcc

Jonathan Cast ha scritto:
On Thu, 2009-02-05 at 01:10 +0100, Manlio Perillo wrote:
Yitzchak Gale ha scritto:
In our case, the Python division first does a quick estimate of the sizes of the two integers, and just returns zero if it sees that there will be underflow on conversion to double. So I made the following rough change to the Haskell:
-- An exact division (/.) :: Integer -> Integer -> Double x /. y | y `div` x > 5*10^323 = 0 | otherwise = fromRational $ toRational x / toRational y
Right, but I would like to see a proper implemented function for exact integer division in GHC.
(%) *is* a proper function for exact integer division. But you'll find plenty of Haskellers to balk at calling anything that returns a Double `proper'.
After a quick search I found the function fromRat in GHC.Float, that implement the division of two big integer returning a floating point number. However it only supports rationals, and with rational there is a performance problem caused by fraction reduction. I would like to do some test with fraction reduction disabled, but right now I don't have the time. Regards Manlio Perillo

On Thu, 5 Feb 2009, Manlio Perillo wrote:
Yitzchak Gale ha scritto:
Ah, OK. Thanks. Now we have a well-defined problem. :)
Good :).
I have used this example to describe how to avoid big integers at all: http://haskell.org/haskellwiki/Integers_too_big_for_floats

Henning Thielemann ha scritto:
On Thu, 5 Feb 2009, Manlio Perillo wrote:
Yitzchak Gale ha scritto:
Ah, OK. Thanks. Now we have a well-defined problem. :)
Good :).
I have used this example to describe how to avoid big integers at all: http://haskell.org/haskellwiki/Integers_too_big_for_floats
Thanks, I didn't know that Ratio constructor was available via GHC.FLoat. I have done a test using the :% operator, instead of the % operator, but, the program run slower! About 29 seconds, against 20 seconds of the previous version (and 10 seconds of the Python version). I was expecting to obtain an execution time of about 6 seconds... Manlio Perillo

Yitzchak Gale schrieb:
Manlio Perillo wrote:
However there is still a *big* problem: it is inefficient.
Here is a Python version of the Chudnovsky algorithm [1] for computing Pi: http://paste.pocoo.org/show/102800/ On my system it takes 10 seconds. Here is an Haskell version: http://paste.pocoo.org/show/102801/ On my system it takes 30 seconds.
Ah, OK. Thanks. Now we have a well-defined problem. :)
You would not need to handle huge interim results, if you compute the factorials incrementally. E.g. (n+1)!^3 = n!^3 * (n+1)^3 etc. 'scanl' is your friend. This should also be faster.

Manlio Perillo wrote:
fac(777) / fac(777) 1.0 Here CPython does not convert the two integers to float before to divide them, but make use of a special algorithm. GHC, instead, returns NaN
I wrote:
No, actually here Haskell shines. Perhaps this GHCi session will illuminate the issue for you: Prelude> let fac n = product [2..n] Prelude> fac 777 `div` fac 777 1 Prelude> fac 777 / fac 777 NaN
No, this is not as complete as it is done in Python. The `/` operator *always* ...return a float... The `//` operator *always* perform a "floor" division. This happens for both integers and floats: If I'm correct, there is no operator/function, in Haskell, that perform an exact division between two integers and return a float: exactDiv :: (Integral a, Real b) => a -> a -> b
exactDiv :: (Real a, Real b, Fractional c) => a -> b -> c x `exactDiv` y = realToFrac x / realToFrac y Python does that same type coercion automatically at runtime. This can sometimes make your code look a bit neater, but other times really gets in the way.
I personally prefer the Python solution, where we have two operators with the same behaviour over all the numbers.
If you really like exactDiv better, you can define an operator to do it for you. But most people don't. Haskell gives you total control, since the types of everything are strictly defined at compile time. This has the extra advantage of giving you a really powerful type checker that can find many of your bugs at compile time. -Yitz

Yitzchak Gale ha scritto:
[...]
exactDiv :: (Real a, Real b, Fractional c) => a -> b -> c x `exactDiv` y = realToFrac x / realToFrac y
Python does that same type coercion automatically at runtime.
No, this is not correct. CPython converts to float only for a "simple" integer (since conversion from machine int to machine float will never overflow). For a "long integer" it uses a special algorithm for computing the division, that avoids overflow if the final result is in the range of a float. http://code.python.org/hg/trunk/file/tip/Objects/longobject.c#l2657
[...]
Thanks Manlio Perillo
participants (9)
-
Henning Thielemann
-
Jonathan Cast
-
Manlio Perillo
-
Max Rabkin
-
Richard O'Keefe
-
Ross Mellgren
-
Simon Peyton-Jones
-
Thomas DuBuisson
-
Yitzchak Gale