
On Nov 14, 2007 12:32 PM, Kurt Hutchinson
As part of a solution I'm working on for Project Euler problem 119, I wanted to create an ordered list of all powers of all positive integers (starting with squares). This is what I came up with:
powers = ( uniq . map fst . iterate next ) ( 1, ( 0, powertable ) )
powertable = map (\ n -> map (\ p -> n ^ p ) [ 2 .. ] ) [ 2 .. ]
next ( _, ( lim, ps ) ) = ( p, ( lim', ps' ) ) where ( c, p ) = minimumBy ( compare `on` snd ) . zip [ 0 .. lim ] . head . transpose $ ps lim' = lim + ( if c == lim then 1 else 0 ) ( front, b : back ) = splitAt c ps ps' = front ++ ( tail b : back )
-- like the unix utility uniq [] = [] uniq [x] = [x] uniq (x:y:ys) | x == y = uniq (y:ys) | otherwise = x : uniq (y:ys)
Basically, think of a grid of numbers, each row is the list of powers for one integer. To find the next power in the sequence, look at the the first number in each row (since that will be the smallest number in the row), up to limit number of rows. The limit is calculated as the number of rows that we've already taken one number from, plus one. This exploits the fact that the row heads are initially ordered. If we use a number from a row, shift that row down to remove the number used.
It does pretty well, but I'd welcome any comments, or even suggestions of a completely different algorithm (for this little exercise, or problem 119). Thanks.
The merging can be done much more simply and efficiently (this is code I wrote when computing squarefree numbers in a blog posthttp://byorgey.wordpress.com/2007/09/01/squarefree-numbers-in-haskell/a few months ago): -- merge two nondecreasing lists. ( # ) :: (Ord a) => [a] -> [a] -> [a] [] # ys = ys xs # [] = xs xs@(x:xt) # ys@(y:yt) | x < y = x : (xt # ys) | x > y = y : (xs # yt) | otherwise = x : (xt # yt) -- merge an infinite list of lists, assuming that each list -- is nondecreasing and the lists occur in order of their first -- element. mergeAll :: (Ord a) => [[a]] -> [a] mergeAll ([] : zs) = mergeAll zs mergeAll (xxs@(x:xs) : yys@(y:ys) : zs) | x < y = x : mergeAll (xs : yys : zs) | otherwise = mergeAll ((xxs # yys) : zs) Then you can simply define powers = 1 : mergeAll powertable I wrote some code to sum the first n powers and print the result, and compiled (using -O2) first with your method, then with mergeAll. Summing the first 7000 powers took ~8s and ~0.1s respectively, so that's a pretty big speedup. Based on seeing how the times scale, I suspect that your code is O(n^2) or something of that sort, whereas mergeAll is essentially linear, although without scrutinizing your code more closely I'm not exactly sure why that would be the case. -Brent