Re: [Haskell-cafe] fast integer base-2 log function?
Hello, If the standard libraries provide such a function, I haven't found it. I must admit that I haven't studied your code in detail. I usually do as follows for integer logarithms, shamelessly stolen from the Haskell report:
-- Integer log base (c.f. Haskell report 14.4):
imLog :: Integer->Integer->Integer imLog b x = if x < b then 0 else let l = 2 * imLog (b*b) x doDiv x l = if x < b then l else doDiv (x`div`b) (l+1) in doDiv (x`div`(b^l)) l
Best regards Thorkil On Monday 11 February 2008 07:15, Uwe Hollerbach wrote:
Hello, haskellers,
Is there a fast integer base-2 log function anywhere in the standard libraries? I wandered through the index, but didn't find anything that looked right. I need something that's more robust than logBase, it needs to handle numbers with a few to many thousands of digits. I found a thread from a couple of years ago that suggested there was no such routine, and that simply doing "length (show n)" might be the best. That seems kind of... less than elegant. I've come up with a routine, shown below, that seems reasonably fast (a few seconds of CPU time for a million-bit number, likely adequate for my purposes), but it seems that something with privileged access to the innards of an Integer ought to be even much faster -- it's just a simple walk along a list (array?) after all. Any pointers? Thanks!
Uwe
powi :: Integer -> Integer -> Integer powi b e | e == 0 = 1 | e < 0 = error "negative exponent in powi" | even e = powi (b*b) (e `quot` 2) | otherwise = b * (powi b (e - 1))
ilog2 :: Integer -> Integer ilog2 n | n < 0 = ilog2 (- n) | n < 2 = 1 | otherwise = up n (1 :: Integer) where up n a = if n < (powi 2 a) then bin (quot a 2) a else up n (2*a) bin lo hi = if (hi - lo) <= 1 then hi else let av = quot (lo + hi) 2 in if n < (powi 2 av) then bin lo av else bin av hi
(This was all properly aligned when I cut'n'pasted; proportional fonts might be messing it up here.) _______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe
Thanks, guys! It looks at first glance as if the code Thorkil posted is similar to mine (grow comparison number in steps of 2 in the exponent, then binary-search to get the exact exponent), while Stefan's version is more similar to the walk-the-list idea I had in mind. I'll play with both of these when I get a chance... Uwe On Feb 10, 2008 10:44 PM, Thorkil Naur <naur@post11.tele.dk> wrote:
Hello,
If the standard libraries provide such a function, I haven't found it. I must admit that I haven't studied your code in detail. I usually do as follows for integer logarithms, shamelessly stolen from the Haskell report:
-- Integer log base (c.f. Haskell report 14.4):
imLog :: Integer->Integer->Integer imLog b x = if x < b then 0 else let l = 2 * imLog (b*b) x doDiv x l = if x < b then l else doDiv (x`div`b) (l+1) in doDiv (x`div`(b^l)) l
Best regards Thorkil
On Monday 11 February 2008 07:15, Uwe Hollerbach wrote:
Hello, haskellers,
Is there a fast integer base-2 log function anywhere in the standard libraries? I wandered through the index, but didn't find anything that looked right. I need something that's more robust than logBase, it needs to handle numbers with a few to many thousands of digits. I found a thread from a couple of years ago that suggested there was no such routine, and that simply doing "length (show n)" might be the best. That seems kind of... less than elegant. I've come up with a routine, shown below, that seems reasonably fast (a few seconds of CPU time for a million-bit number, likely adequate for my purposes), but it seems that something with privileged access to the innards of an Integer ought to be even much faster -- it's just a simple walk along a list (array?) after all. Any pointers? Thanks!
Uwe
powi :: Integer -> Integer -> Integer powi b e | e == 0 = 1 | e < 0 = error "negative exponent in powi" | even e = powi (b*b) (e `quot` 2) | otherwise = b * (powi b (e - 1))
ilog2 :: Integer -> Integer ilog2 n | n < 0 = ilog2 (- n) | n < 2 = 1 | otherwise = up n (1 :: Integer) where up n a = if n < (powi 2 a) then bin (quot a 2) a else up n (2*a) bin lo hi = if (hi - lo) <= 1 then hi else let av = quot (lo + hi) 2 in if n < (powi 2 av) then bin lo av else bin av hi
(This was all properly aligned when I cut'n'pasted; proportional fonts might be messing it up here.) _______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe
Hi, all, a few days ago I had asked about fast integer log-2 routines, and got a couple of answers. I've now had a chance to play with the routines, and here's what I found. Initially, Thorkil's routine taken from the Haskell report was about 30% or so faster than mine. When I replaced the calls to my routine "powi" with native calls to ^, my routine became about 10% faster than Thorkil's routine. (I undoubtedly had some reason for using my own version of powi, but I have no idea anymore what that reason was... :-/ ) I initially thought that that speed difference might be due to the fact that my routine had base 2 hard-wired, whereas his routine is for general bases, but that seems not to be the case: when I modified my version to also do general bases, it stayed pretty much the same. I didn't do enough statistics-gathering to really be absolutely positively certain that my routine is indeed 10.000% faster, but there did seem to be a slight edge in speed there. Here's the latest version, in case anyone's interested. I had previously had effectively a bit-length version; since I made it general base-b I changed it to a log function.
ilogb :: Integer -> Integer -> Integer ilogb b n | n < 0 = ilogb b (- n) | n < b = 0 | otherwise = (up b n 1) - 1 where up b n a = if n < (b ^ a) then bin b (quot a 2) a else up b n (2*a) bin b lo hi = if (hi - lo) <= 1 then hi else let av = quot (lo + hi) 2 in if n < (b ^ av) then bin b lo av else bin b av hi
Stefan's routine is, as expected, much much faster still: I tested the first two routines on numbers with 5 million or so bits and they took ~20 seconds of CPU time, whereas I tested Stefan's routine with numbers with 50 million bits, and it took ~11 seconds of CPU time. The limitation of Stefan's routine is of course that it's limited to base 2 -- it is truly a bit-length routine -- and I guess another potential limitation is that it uses GHC extensions, not pure Haskell (at least I had to compile it with -fglasgow-exts). But it's the speed king if those limitations aren't a problem! Uwe
On Thu, Feb 14, 2008 at 08:23:29PM -0800, Uwe Hollerbach wrote:
Stefan's routine is, as expected, much much faster still: I tested the first two routines on numbers with 5 million or so bits and they took ~20 seconds of CPU time, whereas I tested Stefan's routine with numbers with 50 million bits, and it took ~11 seconds of CPU time. The limitation of Stefan's routine is of course that it's limited to base 2 -- it is truly a bit-length routine -- and I guess another potential limitation is that it uses GHC extensions, not pure Haskell (at least I had to compile it with -fglasgow-exts). But it's the speed king if those limitations aren't a problem!
At the risk of stating the obvious, it also hard-codes 32 bit words and certain configurable (but nobody bothers) details of the GMP data representation (big endian, 0 nails). Stefan
On Thu, Feb 14, 2008 at 8:23 PM, Uwe Hollerbach <uhollerbach@gmail.com> wrote:
Stefan's routine is, as expected, much much faster still: I tested the first two routines on numbers with 5 million or so bits and they took ~20 seconds of CPU time, whereas I tested Stefan's routine with numbers with 50 million bits, and it took ~11 seconds of CPU time.
This seems wrong to me; that routine should take a small constant amount of time. I suspect you are measuring the time to construct the 50-million bit numbers as well. If you constructed a single number and called this routine on it several times I am sure you would get far different results, with the first routines taking ~7-11s each and Stefan's GHC/GMP-magic taking almost nothing. -- ryan
Yes, I suspect you are right. I didn't look into that in much detail, although I did try exchanging "(2 ^ 50000000)" with "(1 `shiftL` 50000000)"; but that didn't make any difference. Uwe On Fri, Feb 15, 2008 at 9:21 AM, Ryan Ingram <ryani.spam@gmail.com> wrote:
On Thu, Feb 14, 2008 at 8:23 PM, Uwe Hollerbach <uhollerbach@gmail.com> wrote:
Stefan's routine is, as expected, much much faster still: I tested the first two routines on numbers with 5 million or so bits and they took ~20 seconds of CPU time, whereas I tested Stefan's routine with numbers with 50 million bits, and it took ~11 seconds of CPU time.
This seems wrong to me; that routine should take a small constant amount of time. I suspect you are measuring the time to construct the 50-million bit numbers as well. If you constructed a single number and called this routine on it several times I am sure you would get far different results, with the first routines taking ~7-11s each and Stefan's GHC/GMP-magic taking almost nothing.
-- ryan
participants (4)
-
Ryan Ingram -
Stefan O'Rear -
Thorkil Naur -
Uwe Hollerbach