CODEWITHSUNDEEP
haskellbooleanmonads
tauros_is_the_best_pokemon
tauros_is_the_best_pokemon

Reputation: 958

List comprehension with IO Bool in Haskell

I'm trying to find probable Mersenne primes by checking Mersenne numbers with the Miller Rabin primality check.

If mersennes is the infinite list of Mersenne numbers, I would like to do something like:

probableMersennePrimes :: IO [Integer]
probableMersennePrimes = [ n | n <- mersennes, primeMillerRabin n ]

Where primeMillerRabin :: Integer -> IO Bool.

As a side question, should the type be IO [Integer] or [IO Integer]?

Upvotes: 1

Views: 1020

Answers (3)

is7s
is7s

Reputation: 3480

All other answers are great. However, I think the simplest approach to this problem is lazy IO and it's worth considering. I also believe that in this special case lazy IO is harmless since there are no system resources involved (file handles etc.)

You just need to redefine a special filterM for lazy IO. 

import System.IO.Unsafe (unsafeInterleaveIO)

filterMIO :: (a -> IO Bool) -> [a] -> IO [a]
filterMIO p = go
  where
    go []     = return []
    go (x:xs) = do
      xs' <- unsafeInterleaveIO (go xs)
      b   <- p x
      return $ if b then (x:xs') else xs'

probableMersennePrimes :: IO [Integer]
probableMersennePrimes = filteMIO primeMillerRabin mersennes

Note that this will just work fine with infinite lists. That's why we've used lazy IO, actually!

Upvotes: 0

rampion
rampion

Reputation: 89073

If you change primeMillerRabin to be of type RandomGen g => Integer -> State g Bool, you can do this with filterM.

probableMersennePrimes :: RandomGen g => g -> [Integer]
probableMersennePrimes = evalState $ filterM primeMillerRabin mersennes

By using evalState, we abandon the final state of the filterM primeMillerRabin mersennes computation, so we cannot be strict in it. But this is a good thing, as the final state will only be available after we reach process the end of the mersennes list, which is infinite, and has no end.

This allows the computation to lazily produce elements of probableMersennePrimes gen.

You can't do it while it returns an IO action, because of how the random number generator works. The IO computation needs to know what the end state is so it can generate another random number after that computation, so it has to loop endlessly, looking for the end of the endless list.

But don't just believe me, try it out:

module SO26307073 where
import Control.Monad.State
import System.Random

-- find how many times a factor divides a number
-- (p^s * d) `factorBy` p == (s,d) iff d `rem` p /= 0
factorBy :: Integral a => a -> a -> (Int,a)
factorBy n p = (length steps - 1, fst $ last steps)
  where steps = takeWhile ((==0) . snd) $ iterate (flip quotRem 2 . fst) (n, 0)

mersennes :: Num a => [a]
mersennes = [ 2^n - 1 | n <- [2..] ]

type RandomRM m = (Integer, Integer) -> m Integer

primeMillerRabinWith :: Monad m => RandomRM m ->  Integer -> m Bool
primeMillerRabinWith randomRM n = do
  let nMinus1 = n-1
      (s,d) = nMinus1 `factorBy` 2
  liftM (all id) . replicateM 10 $ do
    a <- randomRM (2, nMinus1)
    let x = (a^d) `mod` n
    let xs = take s $ iterate (\x -> (x^2) `mod` n) x
    return $ x == 1 || any (== nMinus1) xs

probableMersennePrimesWith :: Monad m => RandomRM m -> m [Integer]
probableMersennePrimesWith randomRM = filterM (primeMillerRabinWith randomRM) mersennes

probableMersennePrimesPure :: RandomGen g => g -> [Integer]
probableMersennePrimesPure = evalState . probableMersennePrimesWith $ state . randomR

probableMersennePrimesIO :: IO [Integer]
probableMersennePrimesIO = probableMersennePrimesWith $ randomRIO

Note that probableMersennePrimesIO and probableMersennePrimesPure just use different ways to draw randoms.

Popping over to ghci, we can see that the pure version works, while the IO version just hangs:

λ import Control.Applicative
λ import System.Random
λ :l SO26307073
[1 of 1] Compiling SO26307073      ( SO26307073.hs, interpreted )
Ok, modules loaded: SO26307073.
λ take 0 . probableMersennePrimesPure <$> newStdGen
Loading package array-0.5.0.0 ... linking ... done.
Loading package deepseq-1.3.0.2 ... linking ... done.
Loading package old-locale-1.0.0.6 ... linking ... done.
Loading package time-1.4.2 ... linking ... done.
Loading package random-1.0.1.1 ... linking ... done.
Loading package transformers-0.4.1.0 ... linking ... done.
Loading package mtl-2.2.1 ... linking ... done.
[]
λ take 5 . probableMersennePrimesPure <$> newStdGen
[3,7,31,127,8191]
λ take 5 <$> probableMersennePrimesIO 
^CInterrupted.

Upvotes: 1

Gabriella Gonzalez
Gabriella Gonzalez

Reputation: 35089

You can stream the infinite list of numbers using ListT from pipes:

import Control.Monad
import Data.Monoid
import Pipes

mersennes :: [Integer]
mersennes = undefined

primeMillerRabin :: Integer -> IO Bool
primeMillerRabin = undefined

probableMersennePrimes :: ListT IO Integer
probableMersennePrimes = do
    n        <- Select (each mersennes)
    continue <- lift (primeMillerRabin n)
    guard continue
    return n

main = runListT $ do
    n <- probableMersennePrimes
    lift (print n)
    mempty

This works even if the list of mersennes is infinite and it will stream results as you compute them.

Upvotes: 4

Related Questions