Dear HCafe-ers,
Yesterday I decided to take a look at the most recent Euler problem,
number 249, and give it a shot. I have a couple of computers at home,
a Dell laptop and a desktop. I compiled this message with ghc -O2
--make ex429.lhs and ran it on each machine. On the Dell I get:
time ./ex429
[650,16900,547924,27396200,746640991,773879749,683631060]
[650,16900,547924,27396200,746640991,773879749,683631060]
136342232
./ex429 8.66s user 0.02s system 99% cpu 8.695 total
When I run this exact same file on the desktop, I get:
time ./ex429
[650,16900,547924,27396200,746640991,773879749,683631060]
[650,16900,547924,27396200,746640991,773879749,683631060]
98792821
./ex429 6.50s user 0.03s system 99% cpu 6.537 total
Which happens to be the right answer. But WHY is the output from the
Dell different?
Machine info is at the bottom of this message.
>{-# LANGUAGE BangPatterns #-}
>{-# LANGUAGE FlexibleContexts #-}
>{-# OPTIONS -O2 -optc-O #-}
>{-
>Sum of squares of unitary divisors
>Problem 429
>A unitary divisor d of a number n is a divisor of n that has the property
>gcd(d, n/d) = 1.
>The unitary divisors of 4! = 24 are 1, 3, 8 and 24.
>The sum of their squares is 12 + 32 + 82 + 242 = 650.
>
>Let S(n) represent the sum of the squares of the unitary divisors of n. Thus
>S(4!)=650.
>
>Find S(100 000 000!) modulo 1 000 000 009.
>
>
>-}
>
>
>import Control.Monad.ST
>import Data.Array.ST
>import Data.Array.IArray as I
>import Data.Array.Unboxed as U
>import Data.Word
>import Data.Ord
>import Data.List
>
>allfactors n = [ i | i <- [1..n] , n `mod` i == 0]
>
>factorial n = product [1..n]
>ud0 n =
> let
> nf = factorial n
> a = allfactors nf
> b = filter (\x -> gcd x (nf `div` x) == 1) a
> in b
>
>ud1 = sum . map (\x -> x*x) . ud0
>ansSlow n = ud1 n `mod` (fromIntegral modulus)
>
>largestExponentInFactorial n p =
> let a = [ n `div` (p^i) | i <- [1..] ]
> b = takeWhile (>0) a
> in sum b
>
>modProduct :: [Int] -> Int
>modProduct = foldl' (\a b -> times a b modulus) 1
>
>pA n = primesA (fromIntegral n)
>
>primesN :: Int -> [Int]
>primesN n = map fromIntegral $ primeS (pA n)
>
>times :: Int -> Int -> Int -> Int
>times x y n =
> let
> x1 = fromIntegral x :: Integer
> y1 = fromIntegral y :: Integer
> n1 = fromIntegral n :: Integer
> result = fromIntegral $! x1 * y1 `mod` n1
> in result
>
>fastPower :: Int -> Int -> Int -> Int
>fastPower x 0 modulus = 1
>fastPower x 1 modulus = x `mod` modulus
>fastPower x n modulus
> | even n = fastPower (times x x modulus) (n `div` 2) modulus
> | otherwise = (times x (fastPower x (n-1) modulus)) modulus
>
>
>foldFun :: Int -> Int -> Int
>foldFun n p =
> let
>a = largestExponentInFactorial n p
>b = fastPower p a modulus
>c = times b b modulus + 1
> in c
>
>ff :: Int -> [Int] -> Int
>ff n = foldl' (\a p -> times a (foldFun n p) modulus) 1
>
>
>ans n =
> let
>ps = primeS $ primesA n
>-- ps = takeWhile (<= n) primes
> in ff n ps
>
>modulus = 19 :: Int
>main = do
> print $ map ans [4..10]
> print $ map ansSlow [4..10]
> print $ ans 1
>
>
>
>{-
>
>intended Usage:
>
>pA = primesA (10^9)
>primes = primeS pA
>isPrime = isPrimE pA
>
>
>
>-}
>
>
>
>sieve :: STUArray s Int Bool -> Int -> Int -> ST s (STUArray s Int Bool)
>sieve !a !m !n
>| n == m= return a
>| otherwise = do
>e <- readArray a (fromIntegral n)
>if e
>then let loop !j
>| j <= m= writeArray a (fromIntegral j) False
> >> loop (j+n)
>| otherwise = sieve a m (n+1)
> in loop (n+n)
>else sieve a m (n+1)
>
>
>primesA :: Int -> UArray Int Bool
>primesA sizeN =
> runSTUArray (do a <- newArray (0,sizeN) True
> :: ST s (STUArray s Int Bool)
> writeArray a 0 False
> writeArray a 1 False
> sieve a sizeN 2)
>
>
>primeS :: (IArray a1 Bool, Ix a) => a1 a Bool -> [a]
>primeS primeArray = map fst $ filter (\x -> snd x) (assocs primeArray)
>
>isPrimE :: (IArray a e, Ix i) => a i e -> i -> e
>isPrimE primeArray n = primeArray I.! n
>
AMD-64 Desktop
uname -a
Linux myth 3.2.0-4-amd64 #1 SMP Debian 3.2.41-2+deb7u2 x86_64 GNU/Linux
ghc --version
The Glorious Glasgow Haskell Compilation System, version 7.6.2
hwinfo --cpu
01: None 00.0: 10103 CPU
[Created at cpu.304]
Unique ID: rdCR.j8NaKXDZtZ6
Hardware Class: cpu
Arch: X86-64
Vendor: "AuthenticAMD"
Model: 21.1.2 "AMD FX(tm)-4100 Quad-Core Processor"
Features:
fpu,vme,de,pse,tsc,msr,pae,mce,cx8,apic,sep,mtrr,pge,mca,cmov,pat,pse36,clflush,mmx,fxsr,sse,sse2,ht,syscall,nx,mmxext,fxsr_opt,pdpe1gb,rdtscp,lm,constant_tsc,rep_good,nopl,nonstop_tsc,extd_apicid,aperfmperf,pni,pclmulqdq,monitor,ssse3,cx16,sse4_1,sse4_2,
Clock: 1400 MHz
BogoMips: 7248.25
Cache: 2048 kb
Units/Proc
