[Haskell-cafe] factorising to prime numbers

[Robert Dockins]

On Feb 9, 2007, at 9:20 AM, Dougal Stanton wrote:

> Hi folks,
>
> I recently read in my copy of Concrete Mathematics the relationship
> between prime factors powers and lcm/gcd functions. So I decided to
> reimplement gcd and lcm the long way, for no other reason than because
> I could.
>
> If you look at the definition of 'powers' you'll note it's  
> infinite. So
> there's no easy way to take the product of this list, if I don't know
> how many items to take from it.
>
> Is there a better way to turn an integer N and a list of primes
> [p1,p2,p3,...] into powers [c1,c2,c3,...] such that
>
> N = product [p1^c1, p2^c2, p3^c3, ...]
>
> If I'm missing something really obvious I'll be very grateful. I can't
> really work out what kind of structure it should be. A map? fold?

If I've understood correctly your list 'powers' will be all zeros  
after a certain point.  Once that happens, you don't need to examine  
that part of the list anymore.  This should at least occur as soon as  
the primes become larger than your number N (and probably sooner.   
sqrt(N) maybe? I forget).  So, you should be able to only examine a  
prefix of the list 'primes'.  The definition you have looks right, in  
that it correctly generates the correct list.  If you want to test  
that its doing the right thing, you can just examine the prefix:

 > test n = product (zipWith (^) (takeWhile (<n) primes) (powers n))

(untested, but I think it would work).

or you can just create the portion of the powers list you need in the  
first place:


 > powersPrefix n = map (f n) (takeWhile (<n) primes)



(remember kids, a decidable problem is a semi-decidable problem where  
we can calculate a stopping condition).


> D.
>
>
> -- Concrete Mathematics
> -- Graham, Knuth & Patashnuk
>
> module Concrete where
>
> import Data.List
>
> -- the sieve of eratosthenes is a fairly simple way
> -- to create a list of prime numbers
> primes =
>     let primes' (n:ns) = n : primes' (filter (\v -> v `mod` n /= 0)  
> ns)
>     in primes' [2..]
>
> -- how many of the prime p are in the unique factorisation
> -- of the integer n?
> f 0 _ = 0
> f n p | n `mod` p == 0 = 1 + f (n `div` p) p
>       | otherwise = 0
>
> powers n = map (f n) primes
>
> --gcd :: Integer -> Integer -> Integer
> --gcd = f . map (uncurry min)
>
> -- 
> Dougal Stanton



Rob Dockins

Speak softly and drive a Sherman tank.
Laugh hard; it's a long way to the bank.
           -- TMBG