[Haskell-cafe] Project Euler Problem 357 in Haskell

[mukesh tiwari]

Hello all
Being a Haskell enthusiastic , first I tried to solve this problem in Haskell but it running for almost 10 minutes on my computer but not getting the answer. A similar C++ program outputs the answer almost instant so could some one please tell me how to improve this Haskell program.

import Control.Monad.ST
import Data.Array.ST
import Data.Array.Unboxed
import Control.Monad

prime :: Int -> UArray Int Bool
prime n = runSTUArray $ do
    arr <- newArray ( 2 , n ) True :: ST s ( STUArray s Int Bool )
    forM_ ( takeWhile ( \x -> x*x <= n ) [ 2 .. n ] ) $ \i -> do
        ai <- readArray arr i
        when ( ai  ) $ forM_ [ i^2 , i^2 + i .. n ] $ \j -> do
            writeArray arr j False

    return arr

pList :: UArray Int Bool
pList = prime $  10 ^ 8
 
divPrime :: Int -> Bool
divPrime n = all ( \d -> if mod n d == 0 then pList ! ( d + div  n  d ) else True )  $  [ 1 .. truncate . sqrt . fromIntegral  $ n ]


main = putStrLn . show . sum  $ [ if and [ pList ! i , divPrime . pred $ i ] then pred  i else 0 | i <- [ 2 .. 10 ^ 8 ] ]


C++ program which outputs the answer almost instant.

#include<cstdio>
#include<iostream>
#include<vector>
#define Lim 100000001
using namespace std;

bool prime [Lim];
vector<int> v ;

void isPrime ()
     {
		for( int i = 2 ; i * i <= Lim ; i++)
		 if ( !prime [i]) for ( int j = i * i ; j <= Lim ; j += i ) prime [j] = 1 ; 

		for( int i = 2 ; i <= Lim ; i++) if ( ! prime[i] ) v.push_back( i ) ;
		//cout<<v.size()<<endl;
		//for(int i=0;i<10;i++) cout<<v[i]<<" ";cout<<endl;

     }

int main()
	{
		isPrime();
		int n = v.size();
		long long sum = 0;
		for(int i = 0 ; i < n ; i ++)
		 {
			int k = v[i]-1;
			bool f = 0;
			for(int i = 1 ; i*i<= k ; i++)
				if ( k % i == 0 && prime[ i + ( k / i ) ] )  { f=1 ; break ; }

			if ( !f ) sum += k;
		 }
		cout<<sum<<endl;
	}

Regards
Mukesh Tiwari

[Lyndon Maydwell]

Could Int be overflowing?

> _______________________________________________
> Haskell-Cafe mailing list
> Haskel...@haskell.org
> http://www.haskell.org/mailman/listinfo/haskell-cafe
>
>

_______________________________________________
Haskell-Cafe mailing list
Haskel...@haskell.org
http://www.haskell.org/mailman/listinfo/haskell-cafe

[mukesh tiwari]

I am  not sure about Int overflow. There is no case of Int overflow in prime , pList and divPrime function however lets assuming Int overflow in main but then still answer should be outputted.

Regards
Mukesh Tiwari

[Ivan Lazar Miljenovic]

May I suggest you try a non-ST solution first (e.g. using Data.IntMap)
first (assuming an auxiliary data-structure is required)?

Also, I'm not sure if the logic in the two versions is the same: I'm
not sure about how you handle the boolean aspect in C++, but you have
a third for-loop there that doesn't seem to correspond to anything in
the Haskell version.

--
Ivan Lazar Miljenovic
Ivan.Mi...@gmail.com
IvanMiljenovic.wordpress.com

[mukesh tiwari]

Logic is same. The idea is generate the primes less than 10^8.  Now from each prime , subtract 1 ( when d is 1 then  d + n / d => n + 1 should be prime ) and check for all the divisor. If all divisor are prime then return True else False

divPrime n = all ( \d -> if mod n d == 0 then pList ! ( d + div  n  d ) else True )  $  [ 1 .. truncate . sqrt . fromIntegral  $ n ]

Corresponding C++ statement

for(int i = 1 ; i*i<= k ; i++) if ( k % i == 0 && prime[ i + ( k / i ) ] )  { f=1 ; break ; }

If al the divPrime returns true then sum this number other wise not. The only difference is after generating the prime in c++ , I  collected all the primes  in vector how ever i don't  think it could be issue for not getting output in Haskell.

On Tue, Nov 8, 2011 at 5:29 PM, Ivan Lazar Miljenovic <ivan.mi...@gmail.com> wrote:

May I suggest you try a non-ST solution first (e.g. using Data.IntMap)
first (assuming an auxiliary data-structure is required)?

Also, I'm not sure if the logic in the two versions is the same: I'm
not sure about how you handle the boolean aspect in C++, but you have
a third for-loop there that doesn't seem to correspond to anything in
the Haskell version.

Which  loop ?
 

[Ivan Lazar Miljenovic]

On 8 November 2011 23:29, mukesh tiwari <mukeshtiw...@gmail.com> wrote:
>> Also, I'm not sure if the logic in the two versions is the same: I'm
>> not sure about how you handle the boolean aspect in C++, but you have
>> a third for-loop there that doesn't seem to correspond to anything in
>> the Haskell version.
>>
> Which  loop ?

for( int i = 2 ; i <= Lim ; i++) if ( ! prime[i] ) v.push_back( i ) ;

--

Ivan Lazar Miljenovic
Ivan.Mi...@gmail.com
IvanMiljenovic.wordpress.com

_______________________________________________

[mukesh tiwari]

In that loop , I  am collecting all the primes in vector how ever I changed the c++ code and  now it resembles to Haskell code. This code  still gives the answer within a second.

#include<cstdio>
#include<iostream>
#include<vector>
#define Lim 100000001
using namespace std;

bool prime [Lim];
vector<int> v ;

void isPrime ()
     {
        for( int i = 2 ; i * i <= Lim ; i++)
         if ( !prime [i]) for ( int j = i * i ; j <= Lim ; j += i ) prime [j] = 1 ;

        //for( int i = 2 ; i <= Lim ; i++) if ( ! prime[i] ) v.push_back( i ) ;

        //cout<<v.size()<<endl;
        //for(int i=0;i<10;i++) cout<<v[i]<<" ";cout<<endl;

     }

int main()
    {
        isPrime();
        int n = v.size();
        long long sum = 0;

        for(int i = 0 ; i < Lim ; i ++)

         if ( ! prime [i])

         {
            int k = i-1;

            bool f = 0;
            for(int i = 1 ; i*i<= k ; i++)
                if ( k % i == 0 && prime[ i + ( k / i ) ] )  { f=1 ; break ; }

            if ( !f ) sum += k;
         }
        cout<<sum<<endl;
    }

Regards
Mukesh Tiwari

[Ryan Yates]

If I compile with optimizations:

ghc --make -O3 primes.hs

I get an answer that is off by one from the C++ program in a few seconds.

[Ryan Yates]

I forgot to add, I'm on 32-bit GHC and the sum will overflow there, so
I changed main:

main = putStrLn . show . sum $ ([ if and [ pList ! i , divPrime .
pred $ i ] then (fromIntegral $ pred i) else 0 | i <- [ 2 .. 10 ^ 8 ]
] :: [Integer])

[mukesh tiwari]

Thank you Ryan .  I never compiled my program with -O3 option before and now i can feel the power of compiler optimisation.
Regards
Mukesh Tiwari

[Silvio Frischknecht]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

On 11/08/2011 02:19 PM, Ryan Yates wrote:
> If I compile with optimizations:
>
> ghc --make -O3 primes.hs
>
> I get an answer that is off by one from the C++ program in a few
> seconds.

nice one. Though i wonder. The problem seems to be that without
optimization sum is not tail-recursive. Is sum meant to not be
tail-recursive?

ghci> sum [1..]

eats up all the memory within seconds while

ghci> foldl' (+) 0 [1..]

does not

So Mukesh if you want your program to run without -Ox you should
probably define your one sum'

import Data.List
sum' = foldl' (+) 0

Silvio
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.11 (GNU/Linux)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org/

iQIcBAEBAgAGBQJOuTSFAAoJEDLsP+zrbatW3W8P/04IPhOqSvAI5Cau+GusInCr
CVu932qNVROMb++NHulEtxKQTS7o/WRqrM5OxKclnOi5Rfi/1Da7ozBKfuLtLSo5
W+bJ8AGYTQAGOoIxbsvhAcCBtA35gColc/56cFUbLexZsd4au6SNRxUe+SdhKDIE
YWPoDW/NU5xJCrW7VtJ8G2qWhkDTOtFWqhp63yG+f6ZbJQpK+j0Z7KVPQ42qUb02
vPddhukb3iqlK990r9/vAeXiiKM9wLA4YQSAgRurmn5R2R++ftq78TWOS9J0H3IF
zspAr19FmEHoHxX3ZiGqyG9thmNwTTz/n2EU4U/Pm070bV2AYKfGeT95XJrp4AkH
Na0wixLJQmN1A22yHbojHkWzykQM8CRlfqiJRZfP/mpYwOSj41/uaZnyyVxD/R/u
BT94XoOEIU/XfGN2l/25O3x6oxWOzhdZ9JVFdeNepdsWHzPFf9oLZIUpyAFRyz7p
0i2Xw5chxpN/kCX0cNOrf0RTo1LqBGcLWmePEL540S3aVKMhfv7Pb/PebWy4nfkI
JMKYiiNmQWs3rYpsX252w8H1hO8R/DpdAF7YvMHAyQ84noTy9B7fbYxN4631Md8G
jG94E7IVOcXx/uQXiMIvJ0P62Bg4Lw5VVLPkOlVPqK36BPcWzsbzXkLCUyIcR0wo
QSbAsBYeUjXtnsUCbhkz
=G+dA
-----END PGP SIGNATURE-----

[Daniel Fischer]

On Tuesday 08 November 2011, 12:21:14, mukesh tiwari wrote:
> Hello all
> Being a Haskell enthusiastic , first I tried to solve this problem in
> Haskell but it running for almost 10 minutes on my computer but not
> getting the answer.

Hmm, finishes in 13.36 seconds here, without any changes.
Of course, it has to be compiled with optimisations, ghc -O2.

> A similar C++ program outputs the answer almost instant

2.85 seconds. g++ -O3.
So, yes, much faster, but not orders of magnitude.

> so could some one please tell me how to improve this Haskell
> program.
>
> import Control.Monad.ST
> import Data.Array.ST
> import Data.Array.Unboxed
> import Control.Monad
>
> prime :: Int -> UArray Int Bool
> prime n = runSTUArray $ do
> arr <- newArray ( 2 , n ) True :: ST s ( STUArray s Int Bool )
> forM_ ( takeWhile ( \x -> x*x <= n ) [ 2 .. n ] ) $ \i -> do
> ai <- readArray arr i
> when ( ai ) $ forM_ [ i^2 , i^2 + i .. n ] $ \j -> do
> writeArray arr j False
>
> return arr

Hmm, would have to look at the core, if the optimiser isn't smart enough to
eliminate the lists, you get considerable overhead from that.

Anyway, readArray/writeArray perform bounds checks, you don't have that in
C++, so if you use unsafeRead and unsafeWrite instead, it'll be faster.
(You're doing the checks in *your* code, no point in repeating it.)

>
> pList :: UArray Int Bool
> pList = prime $ 10 ^ 8
>
> divPrime :: Int -> Bool
> divPrime n = all ( \d -> if mod n d == 0 then pList ! ( d + div n d )
> else True ) $ [ 1 .. truncate . sqrt . fromIntegral $ n ]

Use rem and quot instead of mod and div.
That doesn't make too much difference here, but it gains a bit.

That allocates a list, if you avoid that and check in a loop, like in C++,
it'll be a bit faster.
And instead of (!), use unsafeAt to omit a redundant bounds-check.

>
>
> main = putStrLn . show . sum $ [ if and [ pList ! i , divPrime . pred $
> i ] then pred i else 0 | i <- [ 2 .. 10 ^ 8 ] ]

Dont use

and [condition1, condition2]

that's more readable and faster if written as

condition1 && condition2

Don't use pred, use (i-1) instead.

And you're gratuitously adding a lot of 0s, filter the list

sum [i | i <- [1 .. 99999999], pList ! (i+1) && divPrime i]

However, you're allocating a lot of list cells here, it will be faster if
you calculate the sum in a loop, like you do in C++.

Eliminating the unnecessary bounds-checks and the intermediate lists, it
runs in 4.3 seconds here, not too bad compared to the C++.

However, use a better algorithm.
As is, for each prime p you do trial division on (p-1). For every (p-1)
satisfying the criterion, you do about sqrt(p-1) divisions, that costs a
lot of time. You can make the factorisation (and hence finding of divisors)
cheap if you slightly modify your sieve.

[Daniel Fischer]

On Tuesday 08 November 2011, 14:54:18, Silvio Frischknecht wrote:
> On 11/08/2011 02:19 PM, Ryan Yates wrote:
> > If I compile with optimizations:
> >
> > ghc --make -O3 primes.hs

So far, -O3 is not different from -O2 (-On gives you -O2 for n > 2).

*Never* compile code you want to use without optimisations.

Compiling without optimisations is strictly for development, when
compilation time matters because of frequent recompilation.
Once things have stabilised, compile them with optimisations only.

> >
> > I get an answer that is off by one from the C++ program in a few
> > seconds.
>
> nice one. Though i wonder. The problem seems to be that without
> optimization sum is not tail-recursive. Is sum meant to not be
> tail-recursive?

Well, it is tail recursive (foldl, basically), but not strict.
So without optimisations you get the worst of boths worlds (tail recursion
means no incremental output, as could be possible with foldr for lazy
number types, not being strict in the accumulator means it builds huge
thunks, like foldr with a function strict in the second argument).

>
> ghci> sum [1..]
>
> eats up all the memory within seconds while
>
> ghci> foldl' (+) 0 [1..]
>
> does not
>
> So Mukesh if you want your program to run without -Ox you should
> probably define your one sum'
>
> import Data.List
> sum' = foldl' (+) 0

That'd help, it would still be dog-slow, though, since optimisation is also
crucial for the sieve.

>
> Silvio